diff --git a/CMakeLists.txt b/CMakeLists.txt index df09b8c49..00f3e9222 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -682,6 +682,7 @@ if(CUDA_FOUND) GpuView gpu_ad gpu_linear_two_phase_material + gpuPvt PROPERTIES LABELS ${gpu_label}) endif() diff --git a/CMakeLists_files.cmake b/CMakeLists_files.cmake index 94c5bb4e4..d8c43fe8b 100644 --- a/CMakeLists_files.cmake +++ b/CMakeLists_files.cmake @@ -457,12 +457,14 @@ if (HAVE_CUDA) ADD_CUDA_OR_HIP_FILE(TEST_SOURCE_FILES tests test_solver_adapter.cpp) ADD_CUDA_OR_HIP_FILE(TEST_SOURCE_FILES tests test_gpu_ad.cu) ADD_CUDA_OR_HIP_FILE(TEST_SOURCE_FILES tests test_gpu_linear_two_phase_material.cu) + ADD_CUDA_OR_HIP_FILE(TEST_SOURCE_FILES tests test_gpuPvt.cu) # for loop providing the flag --expt-relaxed-constexpr to fix some cuda issues with constexpr if(NOT CONVERT_CUDA_TO_HIP) set(CU_FILES_NEEDING_RELAXED_CONSTEXPR tests/gpuistl/test_gpu_ad.cu tests/gpuistl/test_gpu_linear_two_phase_material.cu + tests/gpuistl/test_gpuPvt.cu ) foreach(file ${CU_FILES_NEEDING_RELAXED_CONSTEXPR}) diff --git a/opm/simulators/linalg/gpuistl/GpuBuffer.cpp b/opm/simulators/linalg/gpuistl/GpuBuffer.cpp index 35f3e002a..b4132ce09 100644 --- a/opm/simulators/linalg/gpuistl/GpuBuffer.cpp +++ b/opm/simulators/linalg/gpuistl/GpuBuffer.cpp @@ -84,22 +84,29 @@ GpuBuffer::resize(size_t newSize) if (newSize < 1) { OPM_THROW(std::invalid_argument, "Setting a GpuBuffer size to a non-positive number is not allowed"); } - // Allocate memory for the new buffer - T* tmpBuffer = nullptr; - OPM_GPU_SAFE_CALL(cudaMalloc(&tmpBuffer, sizeof(T) * newSize)); - // Move the data from the old to the new buffer with truncation - size_t sizeOfMove = std::min({m_numberOfElements, newSize}); - OPM_GPU_SAFE_CALL(cudaMemcpy(tmpBuffer, - m_dataOnDevice, - sizeOfMove * sizeof(T), - cudaMemcpyDeviceToDevice)); + if (m_numberOfElements == 0) { + // We have no data, so we can just allocate new memory + OPM_GPU_SAFE_CALL(cudaMalloc(&m_dataOnDevice, sizeof(T) * newSize)); + } + else { + // Allocate memory for temporary buffer + T* tmpBuffer = nullptr; + OPM_GPU_SAFE_CALL(cudaMalloc(&tmpBuffer, sizeof(T) * m_numberOfElements)); - // free the old buffer - OPM_GPU_SAFE_CALL(cudaFree(m_dataOnDevice)); + // Move the data from the old to the new buffer with truncation + size_t sizeOfMove = std::min({m_numberOfElements, newSize}); + OPM_GPU_SAFE_CALL(cudaMemcpy(tmpBuffer, + m_dataOnDevice, + sizeOfMove * sizeof(T), + cudaMemcpyDeviceToDevice)); - // swap the buffers - m_dataOnDevice = tmpBuffer; + // free the old buffer + OPM_GPU_SAFE_CALL(cudaFree(m_dataOnDevice)); + + // swap the buffers + m_dataOnDevice = tmpBuffer; + } // update size m_numberOfElements = newSize; diff --git a/opm/simulators/linalg/gpuistl/GpuBuffer.hpp b/opm/simulators/linalg/gpuistl/GpuBuffer.hpp index c08a14279..7272d3c28 100644 --- a/opm/simulators/linalg/gpuistl/GpuBuffer.hpp +++ b/opm/simulators/linalg/gpuistl/GpuBuffer.hpp @@ -81,9 +81,9 @@ public: explicit GpuBuffer(const std::vector& data); /** - * @brief Default constructor that will initialize cublas and allocate 0 bytes of memory + * @brief Default constructor not allocating any memory */ - explicit GpuBuffer(); + GpuBuffer() = default; /** * @brief GpuBuffer allocates new GPU memory of size numberOfElements * sizeof(T) @@ -265,7 +265,7 @@ public: private: T* m_dataOnDevice = nullptr; - size_t m_numberOfElements; + size_t m_numberOfElements = 0; void assertSameSize(const GpuBuffer& other) const; void assertSameSize(size_t size) const; diff --git a/tests/gpuistl/test_gpuPvt.cu b/tests/gpuistl/test_gpuPvt.cu new file mode 100644 index 000000000..e5b703fb2 --- /dev/null +++ b/tests/gpuistl/test_gpuPvt.cu @@ -0,0 +1,269 @@ +#include + +#define BOOST_TEST_MODULE TestGpuPvt + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include +#include + +#include +#include +#include +#include + +using Evaluation = Opm::DenseAd::Evaluation; +using GpuB = const Opm::gpuistl::GpuBuffer; +using GpuV = Opm::gpuistl::GpuView; + +using GpuTab = Opm::UniformTabulated2DFunction; + +using GpuBufCo2Tables = Opm::CO2Tables; +using GpuViewCO2Tables = Opm::CO2Tables; +using GpuCO2 = Opm::CO2; + +using HuDuan = Opm::SimpleHuDuanH2O; +using BrineDyn = Opm::BrineDynamic; + +using CpuBrine_CO2 = Opm::BinaryCoeff::Brine_CO2>; +using GpuBrine_CO2 = Opm::BinaryCoeff::Brine_CO2; + +using CpuCo2Pvt = Opm::Co2GasPvt; +using GpuBufCo2Pvt = Opm::Co2GasPvt; +using GpuViewCo2Pvt = Opm::Co2GasPvt; + +namespace { + +/* + This file contains unit tests for Pvt objects and function on the GPU, with additional helper classes +*/ + +const double ABS_TOL = 1e-6; + +struct Fixture { + Fixture(){ + temp = Evaluation(290.5); + pressure = Evaluation(200000.0); + + gpuComputedResultOnCpu = 0.0; + + // move pvt evaluations to gpu + OPM_GPU_SAFE_CALL(cudaMalloc(&gpuTemp, sizeof(Evaluation))); + OPM_GPU_SAFE_CALL(cudaMemcpy(gpuTemp, &temp, sizeof(Evaluation), cudaMemcpyHostToDevice)); + OPM_GPU_SAFE_CALL(cudaMalloc(&gpuPressure, sizeof(Evaluation))); + OPM_GPU_SAFE_CALL(cudaMemcpy(gpuPressure, &pressure, sizeof(Evaluation), cudaMemcpyHostToDevice)); + } + ~Fixture(){ + OPM_GPU_SAFE_CALL(cudaFree(gpuTemp)); + OPM_GPU_SAFE_CALL(cudaFree(gpuPressure)); + } + + Evaluation temp; // [K] + Evaluation pressure; // [Pa] + Evaluation* gpuTemp; // [K] + Evaluation* gpuPressure; // [Pa] + + double gpuComputedResultOnCpu; + + Opm::CO2Tables> co2Tables; +}; + +// Kernel to evaluate a 2D function on the GPU +__global__ void gpuEvaluateUniformTabulated2DFunction(GpuTab gpuTab, Evaluation* inputX, Evaluation* inputY, double* result) { + *result = gpuTab.eval(*inputX, *inputY, true).value(); +} + +// Kernel using a CO2 object on the GPU +__global__ void gpuCO2GasViscosity(GpuViewCO2Tables gpuViewCo2Tables, Evaluation* temp, Evaluation* pressure, double* result) { + *result = GpuCO2::gasViscosity(gpuViewCo2Tables, *temp, *pressure, true).value(); +} + +// Kernel using a SimpleHuDuanH20 object on a GPU +__global__ void huDuanLiquidDensity(Evaluation* temp, Evaluation* pressure, double* result) { + *result = HuDuan::liquidDensity(*temp, *pressure, true).value(); +} + +// Kernel using a BrineDynamic object on a GPU +__global__ void brineDynamicLiquidEnthalpy(Evaluation* temp, Evaluation* pressure, Evaluation* salinity, double* result) { + *result = BrineDyn::liquidEnthalpy(*temp, *pressure, *salinity).value(); +} + +// Kernel using a Brine_CO2 object on a GPU +__global__ void brineCO2GasDiffCoeff(GpuViewCO2Tables co2tables, Evaluation* temp, Evaluation* pressure, double* result) { + *result = GpuBrine_CO2::gasDiffCoeff(co2tables, *temp, *pressure, true).value(); +} + +// Kernel using a Co2GasPvt object on a GPU +__global__ void co2GasPvtInternalEnergy(GpuViewCo2Pvt gpuViewCo2Pvt, Evaluation* temp, Evaluation* pressure, double* result) { + *result = gpuViewCo2Pvt.internalEnergy(1, *temp, *pressure, Evaluation(0.4), Evaluation(0.0)).value(); +} + +// Helper function to launch a kernel and retrieve the result on the CPU to reduce code duplicatoin +template +double launchKernelAndRetrieveResult(KernelFunc kernel, Args... args) { + double* resultOnGpu; + double gpuComputedResultOnCpu; + + // Allocate memory for the result on the GPU + OPM_GPU_SAFE_CALL(cudaMalloc(&resultOnGpu, sizeof(double))); + + // Launch the kernel + kernel<<<1, 1>>>(args..., resultOnGpu); + + // Check for any errors in kernel launch + OPM_GPU_SAFE_CALL(cudaPeekAtLastError()); + OPM_GPU_SAFE_CALL(cudaDeviceSynchronize()); + + // Retrieve the result from the GPU to the CPU + OPM_GPU_SAFE_CALL(cudaMemcpy(&gpuComputedResultOnCpu, resultOnGpu, sizeof(double), cudaMemcpyDeviceToHost)); + + // Free allocated GPU memory + OPM_GPU_SAFE_CALL(cudaFree(resultOnGpu)); + + return gpuComputedResultOnCpu; +} + +bool compareSignificantDigits(double a, double b, int significantDigits) { + // Handle the case where both values are exactly equal + if (a == b) { + return true; + } + + // Calculate the relative error + double relativeError = std::abs(a - b) / std::max(std::abs(a), std::abs(b)); + + // Compute the number of matching digits + double digitsMatched = -std::log10(relativeError); + + // Return true if they match the required number of significant digits + return digitsMatched >= significantDigits; +} + +} // END EMPTY NAMESPACE + +// Test case for evaluating a tabulated 2D function on both CPU and GPU +BOOST_FIXTURE_TEST_CASE(TestEvaluateUniformTabulated2DFunctionOnGpu, Fixture) { + // Example tabulated data (2D) + std::vector> tabData = {{1.0, 2.0}, {3.0, 4.0}, {5.0, 6.0}}; + + // CPU-side function definition + Opm::UniformTabulated2DFunction cpuTab(1.0, 6.0, 3, 1.0, 6.0, 2, tabData); + + // Move data to GPU buffer and create a view for GPU operations + Opm::UniformTabulated2DFunction gpuBufTab = Opm::gpuistl::move_to_gpu(cpuTab); + GpuTab gpuViewTab = Opm::gpuistl::make_view(gpuBufTab); + + // Evaluation points on the CPU + Evaluation a(2.3); + Evaluation b(4.5); + + // Allocate GPU memory for the Evaluation inputs + Evaluation* gpuA = nullptr; + Evaluation* gpuB = nullptr; + OPM_GPU_SAFE_CALL(cudaMalloc(&gpuA, sizeof(Evaluation))); + OPM_GPU_SAFE_CALL(cudaMemcpy(gpuA, &a, sizeof(Evaluation), cudaMemcpyHostToDevice)); + OPM_GPU_SAFE_CALL(cudaMalloc(&gpuB, sizeof(Evaluation))); + OPM_GPU_SAFE_CALL(cudaMemcpy(gpuB, &b, sizeof(Evaluation), cudaMemcpyHostToDevice)); + + gpuComputedResultOnCpu = launchKernelAndRetrieveResult(gpuEvaluateUniformTabulated2DFunction, gpuViewTab, gpuA, gpuB); + + // Free allocated GPU memory + OPM_GPU_SAFE_CALL(cudaFree(gpuA)); + OPM_GPU_SAFE_CALL(cudaFree(gpuB)); + + // Verify that the CPU and GPU results match within a reasonable tolerance + const double cpuComputedResult = cpuTab.eval(a, b, true).value(); + BOOST_CHECK(std::fabs(gpuComputedResultOnCpu - cpuComputedResult) < ABS_TOL); +} + +// Test case evaluating CO2 pvt properties on CPU and GPU +BOOST_FIXTURE_TEST_CASE(TestUseCO2OnGpu, Fixture) { + + // use the CO2 tables to aquire the viscosity at 290[K] and 2e5[Pa] + double viscosityReference = Opm::CO2>>::gasViscosity(co2Tables, temp, pressure, true).value(); + + GpuBufCo2Tables gpuBufCo2Table = Opm::gpuistl::move_to_gpu, GpuB>(co2Tables); + GpuViewCO2Tables gpuViewCo2Table = Opm::gpuistl::make_view(gpuBufCo2Table); + + gpuComputedResultOnCpu = launchKernelAndRetrieveResult(gpuCO2GasViscosity, gpuViewCo2Table, gpuTemp, gpuPressure); + + // Verify that the CPU and GPU results match within a reasonable tolerance + BOOST_CHECK(std::fabs(gpuComputedResultOnCpu - viscosityReference) < ABS_TOL); +} + +// Test case evaluating pvt values for SimpleHuDuanH20 on a GPU and CPU +BOOST_FIXTURE_TEST_CASE(TestUseH2OOnGpu, Fixture) { + + // use the CO2 tables to aquire the densityReference at 290[K] and 2e5[Pa] + double densityReference = HuDuan::liquidDensity(temp, pressure, true).value(); + + gpuComputedResultOnCpu = launchKernelAndRetrieveResult(huDuanLiquidDensity, gpuTemp, gpuPressure); + + // Verify that the CPU and GPU results match within a reasonable tolerance + BOOST_CHECK(std::fabs(gpuComputedResultOnCpu - densityReference) < ABS_TOL); +} + +// Test case evaluating pvt values for BrineDynamic on a GPU and CPU +BOOST_FIXTURE_TEST_CASE(TestUseBrineDynamicOnGpu, Fixture) { + Evaluation salinity(0.1); // [g/Kg] + + // use the CO2 tables to aquire the enthalpyReference at 290[K] and 2e5[Pa] + double enthalpyReference = BrineDyn::liquidEnthalpy(temp, pressure, salinity).value(); + + // Allocate GPU memory for the Evaluation inputs + Evaluation* gpuSalinity = nullptr; + OPM_GPU_SAFE_CALL(cudaMalloc(&gpuSalinity, sizeof(Evaluation))); + OPM_GPU_SAFE_CALL(cudaMemcpy(gpuSalinity, &salinity, sizeof(Evaluation), cudaMemcpyHostToDevice)); + + gpuComputedResultOnCpu = launchKernelAndRetrieveResult(brineDynamicLiquidEnthalpy, gpuTemp, gpuPressure, gpuSalinity); + + // Verify that the CPU and GPU results match within a reasonable tolerance + BOOST_CHECK(std::fabs(gpuComputedResultOnCpu - enthalpyReference) < ABS_TOL); +} + +// Test case evaluating pvt values for BrineDynamic on a GPU and CPU +BOOST_FIXTURE_TEST_CASE(TestBrine_CO2OnGPU, Fixture) { + + // use the CO2 tables to aquire the gasDiffCoeffReference at 290[K] and 2e5[Pa] + double gasDiffCoeffReference = CpuBrine_CO2::gasDiffCoeff(co2Tables, temp, pressure, true).value(); + + // use the CO2 tables to aquire the viscosity at 290[K] and 2e5[Pa] + double viscosity = Opm::CO2>>::gasViscosity(co2Tables, temp, pressure, true).value(); + + GpuBufCo2Tables gpuBufCo2Table = Opm::gpuistl::move_to_gpu, GpuB>(co2Tables); + GpuViewCO2Tables gpuViewCo2Table = Opm::gpuistl::make_view(gpuBufCo2Table); + + gpuComputedResultOnCpu = launchKernelAndRetrieveResult(brineCO2GasDiffCoeff, gpuViewCo2Table, gpuTemp, gpuPressure); + + // Verify that the CPU and GPU results match within a reasonable tolerance + BOOST_CHECK(std::fabs(gpuComputedResultOnCpu - gasDiffCoeffReference) < ABS_TOL); +} + +// Test case evaluating pvt values for BrineDynamic on a GPU and CPU +BOOST_FIXTURE_TEST_CASE(TestCo2GasPvt, Fixture) { + std::vector salinities = {0.2, 0.3, 0.4}; + + CpuCo2Pvt cpuCo2Pvt(salinities); + double internalEnergyReference = cpuCo2Pvt.internalEnergy(1, temp, pressure, Evaluation(0.4), Evaluation(0.0)).value(); + + const GpuBufCo2Pvt gpuBufCo2Pvt = Opm::gpuistl::move_to_gpu(cpuCo2Pvt); + const auto brineReferenceDensityCPUCopy = gpuBufCo2Pvt.getBrineReferenceDensity().asStdVector(); + const GpuViewCo2Pvt gpuViewCo2Pvt = Opm::gpuistl::make_view(gpuBufCo2Pvt); + + gpuComputedResultOnCpu = launchKernelAndRetrieveResult(co2GasPvtInternalEnergy, gpuViewCo2Pvt, gpuTemp, gpuPressure); + + // Verify that the CPU and GPU results match within a reasonable tolerance + BOOST_CHECK(compareSignificantDigits(gpuComputedResultOnCpu, internalEnergyReference, 6)); +}