cuda_vector.cc 5.5 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185
  1. // Ceres Solver - A fast non-linear least squares minimizer
  2. // Copyright 2023 Google Inc. All rights reserved.
  3. // http://ceres-solver.org/
  4. //
  5. // Redistribution and use in source and binary forms, with or without
  6. // modification, are permitted provided that the following conditions are met:
  7. //
  8. // * Redistributions of source code must retain the above copyright notice,
  9. // this list of conditions and the following disclaimer.
  10. // * Redistributions in binary form must reproduce the above copyright notice,
  11. // this list of conditions and the following disclaimer in the documentation
  12. // and/or other materials provided with the distribution.
  13. // * Neither the name of Google Inc. nor the names of its contributors may be
  14. // used to endorse or promote products derived from this software without
  15. // specific prior written permission.
  16. //
  17. // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
  18. // AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
  19. // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
  20. // ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE
  21. // LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
  22. // CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
  23. // SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
  24. // INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
  25. // CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
  26. // ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
  27. // POSSIBILITY OF SUCH DAMAGE.
  28. //
  29. // Author: joydeepb@cs.utexas.edu (Joydeep Biswas)
  30. //
  31. // A simple CUDA vector class.
  32. // This include must come before any #ifndef check on Ceres compile options.
  33. // clang-format off
  34. #include "ceres/internal/config.h"
  35. // clang-format on
  36. #include <math.h>
  37. #include "ceres/context_impl.h"
  38. #include "ceres/internal/export.h"
  39. #include "ceres/types.h"
  40. #ifndef CERES_NO_CUDA
  41. #include "ceres/cuda_buffer.h"
  42. #include "ceres/cuda_kernels_vector_ops.h"
  43. #include "ceres/cuda_vector.h"
  44. #include "cublas_v2.h"
  45. namespace ceres::internal {
  46. CudaVector::CudaVector(ContextImpl* context, int size)
  47. : context_(context), data_(context, size) {
  48. DCHECK_NE(context, nullptr);
  49. DCHECK(context->IsCudaInitialized());
  50. Resize(size);
  51. }
  52. CudaVector::CudaVector(CudaVector&& other)
  53. : num_rows_(other.num_rows_),
  54. context_(other.context_),
  55. data_(std::move(other.data_)),
  56. descr_(other.descr_) {
  57. other.num_rows_ = 0;
  58. other.descr_ = nullptr;
  59. }
  60. CudaVector& CudaVector::operator=(const CudaVector& other) {
  61. if (this != &other) {
  62. Resize(other.num_rows());
  63. data_.CopyFromGPUArray(other.data_.data(), num_rows_);
  64. }
  65. return *this;
  66. }
  67. void CudaVector::DestroyDescriptor() {
  68. if (descr_ != nullptr) {
  69. CHECK_EQ(cusparseDestroyDnVec(descr_), CUSPARSE_STATUS_SUCCESS);
  70. descr_ = nullptr;
  71. }
  72. }
  73. CudaVector::~CudaVector() { DestroyDescriptor(); }
  74. void CudaVector::Resize(int size) {
  75. data_.Reserve(size);
  76. num_rows_ = size;
  77. DestroyDescriptor();
  78. CHECK_EQ(cusparseCreateDnVec(&descr_, num_rows_, data_.data(), CUDA_R_64F),
  79. CUSPARSE_STATUS_SUCCESS);
  80. }
  81. double CudaVector::Dot(const CudaVector& x) const {
  82. double result = 0;
  83. CHECK_EQ(cublasDdot(context_->cublas_handle_,
  84. num_rows_,
  85. data_.data(),
  86. 1,
  87. x.data(),
  88. 1,
  89. &result),
  90. CUBLAS_STATUS_SUCCESS)
  91. << "CuBLAS cublasDdot failed.";
  92. return result;
  93. }
  94. double CudaVector::Norm() const {
  95. double result = 0;
  96. CHECK_EQ(cublasDnrm2(
  97. context_->cublas_handle_, num_rows_, data_.data(), 1, &result),
  98. CUBLAS_STATUS_SUCCESS)
  99. << "CuBLAS cublasDnrm2 failed.";
  100. return result;
  101. }
  102. void CudaVector::CopyFromCpu(const double* x) {
  103. data_.CopyFromCpu(x, num_rows_);
  104. }
  105. void CudaVector::CopyFromCpu(const Vector& x) {
  106. if (x.rows() != num_rows_) {
  107. Resize(x.rows());
  108. }
  109. CopyFromCpu(x.data());
  110. }
  111. void CudaVector::CopyTo(Vector* x) const {
  112. CHECK(x != nullptr);
  113. x->resize(num_rows_);
  114. data_.CopyToCpu(x->data(), num_rows_);
  115. }
  116. void CudaVector::CopyTo(double* x) const {
  117. CHECK(x != nullptr);
  118. data_.CopyToCpu(x, num_rows_);
  119. }
  120. void CudaVector::SetZero() {
  121. // Allow empty vector to be zeroed
  122. if (num_rows_ == 0) return;
  123. CHECK(data_.data() != nullptr);
  124. CudaSetZeroFP64(data_.data(), num_rows_, context_->DefaultStream());
  125. }
  126. void CudaVector::Axpby(double a, const CudaVector& x, double b) {
  127. if (&x == this) {
  128. Scale(a + b);
  129. return;
  130. }
  131. CHECK_EQ(num_rows_, x.num_rows_);
  132. if (b != 1.0) {
  133. // First scale y by b.
  134. CHECK_EQ(
  135. cublasDscal(context_->cublas_handle_, num_rows_, &b, data_.data(), 1),
  136. CUBLAS_STATUS_SUCCESS)
  137. << "CuBLAS cublasDscal failed.";
  138. }
  139. // Then add a * x to y.
  140. CHECK_EQ(cublasDaxpy(context_->cublas_handle_,
  141. num_rows_,
  142. &a,
  143. x.data(),
  144. 1,
  145. data_.data(),
  146. 1),
  147. CUBLAS_STATUS_SUCCESS)
  148. << "CuBLAS cublasDaxpy failed.";
  149. }
  150. void CudaVector::DtDxpy(const CudaVector& D, const CudaVector& x) {
  151. CudaDtDxpy(
  152. data_.data(), D.data(), x.data(), num_rows_, context_->DefaultStream());
  153. }
  154. void CudaVector::Scale(double s) {
  155. CHECK_EQ(
  156. cublasDscal(context_->cublas_handle_, num_rows_, &s, data_.data(), 1),
  157. CUBLAS_STATUS_SUCCESS)
  158. << "CuBLAS cublasDscal failed.";
  159. }
  160. } // namespace ceres::internal
  161. #endif // CERES_NO_CUDA