|
| 1 | +// Copyright John Maddock 2016. |
| 2 | +// Copyright Matt Borland 2024 - 2026. |
| 3 | +// Use, modification and distribution are subject to the |
| 4 | +// Boost Software License, Version 1.0. (See accompanying file |
| 5 | +// LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt) |
| 6 | + |
| 7 | +#ifdef BOOST_MATH_ENABLE_CUDA |
| 8 | +# undef BOOST_MATH_ENABLE_CUDA |
| 9 | +#endif // BOOST_MATH_ENABLE_CUDA |
| 10 | + |
| 11 | +// Purposefully pull in headers that caused errors in the linked issue |
| 12 | +#include <boost/math/special_functions.hpp> |
| 13 | +#include <boost/math/tools/config.hpp> |
| 14 | +#include <boost/math/tools/numeric_limits.hpp> |
| 15 | +#include <boost/math/constants/constants.hpp> |
| 16 | + |
| 17 | +#include <iostream> |
| 18 | +#include <iomanip> |
| 19 | +#include <vector> |
| 20 | +#include <cmath> |
| 21 | +#include <random> |
| 22 | +#include "cuda_managed_ptr.hpp" |
| 23 | +#include "stopwatch.hpp" |
| 24 | + |
| 25 | +#ifdef BOOST_MATH_ENABLE_CUDA |
| 26 | +# error "We should not be enabling this ourselves" |
| 27 | +#endif // BOOST_MATH_ENABLE_CUDA |
| 28 | + |
| 29 | +// For the CUDA runtime routines (prefixed with "cuda_") |
| 30 | +#include <cuda_runtime.h> |
| 31 | + |
| 32 | +typedef double float_type; |
| 33 | + |
| 34 | +/** |
| 35 | + * CUDA Kernel Device code |
| 36 | + * |
| 37 | + */ |
| 38 | +__global__ void cuda_test(const float_type *in, float_type *out, int numElements) |
| 39 | +{ |
| 40 | + using std::cos; |
| 41 | + const int i = blockDim.x * blockIdx.x + threadIdx.x; |
| 42 | + |
| 43 | + if (i < numElements) |
| 44 | + { |
| 45 | + out[i] = cos(in[i]); |
| 46 | + } |
| 47 | +} |
| 48 | + |
| 49 | +/** |
| 50 | + * Host main routine |
| 51 | + */ |
| 52 | +int main() |
| 53 | +{ |
| 54 | + // Error code to check return values for CUDA calls |
| 55 | + cudaError_t err = cudaSuccess; |
| 56 | + |
| 57 | + // Print the vector length to be used, and compute its size |
| 58 | + int numElements = 50000; |
| 59 | + std::cout << "[Vector operation on " << numElements << " elements]" << std::endl; |
| 60 | + |
| 61 | + // Allocate the managed input vector A |
| 62 | + cuda_managed_ptr<float_type> input_vector(numElements); |
| 63 | + |
| 64 | + // Allocate the managed output vector C |
| 65 | + cuda_managed_ptr<float_type> output_vector(numElements); |
| 66 | + |
| 67 | + // Initialize the input vectors |
| 68 | + // Check some of our numeric_limits for viability |
| 69 | + std::mt19937_64 rng {42}; |
| 70 | + std::uniform_real_distribution<float_type> dist(0, boost::math::constants::pi<float_type>()); |
| 71 | + static_assert(boost::math::numeric_limits<float_type>::is_specialized, "Should be since it's a double"); |
| 72 | + static_assert(boost::math::numeric_limits<float_type>::is_signed, "Should be since it's a double"); |
| 73 | + |
| 74 | + for (int i = 0; i < numElements; ++i) |
| 75 | + { |
| 76 | + input_vector[i] = dist(rng); |
| 77 | + } |
| 78 | + |
| 79 | + // Launch the Vector Add CUDA Kernel |
| 80 | + int threadsPerBlock = 256; |
| 81 | + int blocksPerGrid =(numElements + threadsPerBlock - 1) / threadsPerBlock; |
| 82 | + std::cout << "CUDA kernel launch with " << blocksPerGrid << " blocks of " << threadsPerBlock << " threads" << std::endl; |
| 83 | + |
| 84 | + watch w; |
| 85 | + |
| 86 | + cuda_test<<<blocksPerGrid, threadsPerBlock>>>(input_vector.get(), output_vector.get(), numElements); |
| 87 | + cudaDeviceSynchronize(); |
| 88 | + |
| 89 | + std::cout << "CUDA kernal done in: " << w.elapsed() << "s" << std::endl; |
| 90 | + |
| 91 | + err = cudaGetLastError(); |
| 92 | + |
| 93 | + if (err != cudaSuccess) |
| 94 | + { |
| 95 | + std::cerr << "Failed to launch vectorAdd kernel (error code " << cudaGetErrorString(err) << ")!" << std::endl; |
| 96 | + return EXIT_FAILURE; |
| 97 | + } |
| 98 | + |
| 99 | + // Verify that the result vector is correct |
| 100 | + std::vector<float_type> results; |
| 101 | + results.reserve(numElements); |
| 102 | + w.reset(); |
| 103 | + for(int i = 0; i < numElements; ++i) |
| 104 | + results.push_back(std::cos(input_vector[i])); |
| 105 | + double t = w.elapsed(); |
| 106 | + // check the results |
| 107 | + for(int i = 0; i < numElements; ++i) |
| 108 | + { |
| 109 | + if (boost::math::epsilon_difference(output_vector[i], results[i]) > 10) |
| 110 | + { |
| 111 | + std::cerr << "Result verification failed at element " << i << "!" << std::endl; |
| 112 | + return EXIT_FAILURE; |
| 113 | + } |
| 114 | + } |
| 115 | + |
| 116 | + std::cout << "Test PASSED, normal calculation time: " << t << "s" << std::endl; |
| 117 | + std::cout << "Done\n"; |
| 118 | + |
| 119 | + return 0; |
| 120 | +} |
0 commit comments