Add support for dense CUDA solvers #3

1. Add CUDA initialization and cleanup management to the
   ContextImpl object. The ContextImpl is now solely responsible
   for managing CUDA-related resources.
2. All CUDA dense solvers now use lazy CUDA initialization
   via the ContextImpl object.

Change-Id: Ief456860c72e462367ee997d389c19e2bff50baf

internal/ceres/context_impl.cc

@@ -30,11 +30,73 @@
 
 #include "ceres/context_impl.h"
 
+#include <string>
+
+#ifndef CERES_NO_CUDA
+#include "cuda_runtime.h"
+#include "cublas_v2.h"
+#include "cusolverDn.h"
+#endif  // CERES_NO_CUDA
+
 namespace ceres {
 namespace internal {
 
 ContextImpl::ContextImpl() = default;
 
+#ifndef CERES_NO_CUDA
+bool ContextImpl::InitCUDA(std::string* message) {
+  if (cuda_initialized_) {
+    return true;
+  }
+  if (cublasCreate(&cublas_handle_) != CUBLAS_STATUS_SUCCESS) {
+    *message = "cuBLAS::cublasCreate failed.";
+    cublas_handle_ = nullptr;
+    return false;
+  }
+  if (cusolverDnCreate(&cusolver_handle_) != CUSOLVER_STATUS_SUCCESS) {
+    *message = "cuSolverDN::cusolverDnCreate failed.";
+    cusolver_handle_ = nullptr;
+    cublasDestroy(cublas_handle_);
+    cublas_handle_ = nullptr;
+    return false;
+  }
+  if (cudaStreamCreateWithFlags(&stream_, cudaStreamNonBlocking) !=
+      cudaSuccess) {
+    *message = "CUDA::cudaStreamCreateWithFlags failed.";
+    cusolverDnDestroy(cusolver_handle_);
+    cublasDestroy(cublas_handle_);
+    cusolver_handle_ = nullptr;
+    cublas_handle_ = nullptr;
+    stream_ = nullptr;
+    return false;
+  }
+  if (cusolverDnSetStream(cusolver_handle_, stream_) !=
+      CUSOLVER_STATUS_SUCCESS ||
+      cublasSetStream(cublas_handle_, stream_) != CUBLAS_STATUS_SUCCESS) {
+    *message =
+        "cuSolverDN::cusolverDnSetStream or cuBLAS::cublasSetStream failed.";
+    cusolverDnDestroy(cusolver_handle_);
+    cublasDestroy(cublas_handle_);
+    cudaStreamDestroy(stream_);
+    cusolver_handle_ = nullptr;
+    cublas_handle_ = nullptr;
+    stream_ = nullptr;
+    return false;
+  }
+  cuda_initialized_ = true;
+  return true;
+}
+#endif  // CERES_NO_CUDA
+
+ContextImpl::~ContextImpl() {
+#ifndef CERES_NO_CUDA
+  if (cuda_initialized_) {
+    cusolverDnDestroy(cusolver_handle_);
+    cublasDestroy(cublas_handle_);
+    cudaStreamDestroy(stream_);
+  }
+#endif  // CERES_NO_CUDA
+}
 void ContextImpl::EnsureMinimumThreads(int num_threads) {
 #ifdef CERES_USE_CXX_THREADS
   thread_pool.Resize(num_threads);

internal/ceres/context_impl.h

@@ -40,6 +40,12 @@
 #include "ceres/internal/disable_warnings.h"
 #include "ceres/internal/export.h"
 
+#ifndef CERES_NO_CUDA
+#include "cuda_runtime.h"
+#include "cublas_v2.h"
+#include "cusolverDn.h"
+#endif  // CERES_NO_CUDA
+
 #ifdef CERES_USE_CXX_THREADS
 #include "ceres/thread_pool.h"
 #endif  // CERES_USE_CXX_THREADS
@@ -50,6 +56,7 @@ namespace internal {
 class CERES_NO_EXPORT ContextImpl : public Context {
  public:
   ContextImpl();
+  ~ContextImpl() override;
   ContextImpl(const ContextImpl&) = delete;
   void operator=(const ContextImpl&) = delete;
 
@@ -62,6 +69,23 @@ class CERES_NO_EXPORT ContextImpl : public Context {
 #ifdef CERES_USE_CXX_THREADS
   ThreadPool thread_pool;
 #endif  // CERES_USE_CXX_THREADS
+
+#ifndef CERES_NO_CUDA
+  // Initializes the cuSolverDN context, creates an asynchronous stream, and
+  // associates the stream with cuSolverDN. Returns true iff initialization was
+  // successful, else it returns false and a human-readable error message is
+  // returned.
+  bool InitCUDA(std::string* message);
+
+  // Handle to the cuSOLVER context.
+  cusolverDnHandle_t cusolver_handle_ = nullptr;
+  // Handle to cuBLAS context.
+  cublasHandle_t cublas_handle_ = nullptr;
+  // CUDA device stream.
+  cudaStream_t stream_ = nullptr;
+  // Indicates whether all the CUDA resources have been initialized.
+  bool cuda_initialized_ = false;
+#endif  // CERES_NO_CUDA
 };
 
 }  // namespace internal

internal/ceres/cuda_dense_cholesky_test.cc

@@ -43,6 +43,8 @@ namespace internal {
 
 TEST(CUDADenseCholesky, InvalidOptionOnCreate) {
   LinearSolver::Options options;
+  ContextImpl context;
+  options.context = &context;
   auto dense_cuda_solver = CUDADenseCholesky::Create(options);
   EXPECT_EQ(dense_cuda_solver, nullptr);
 }
@@ -56,6 +58,8 @@ TEST(CUDADenseCholesky, Cholesky4x4Matrix) {
         0,   0,   0, 1;
   const Eigen::Vector4d b = Eigen::Vector4d::Ones();
   LinearSolver::Options options;
+  ContextImpl context;
+  options.context = &context;
   options.dense_linear_algebra_library_type = CUDA;
   auto dense_cuda_solver = CUDADenseCholesky::Create(options);
   ASSERT_NE(dense_cuda_solver, nullptr);
@@ -80,6 +84,8 @@ TEST(CUDADenseCholesky, SingularMatrix) {
         0, 0, 0;
   const Eigen::Vector3d b = Eigen::Vector3d::Ones();
   LinearSolver::Options options;
+  ContextImpl context;
+  options.context = &context;
   options.dense_linear_algebra_library_type = CUDA;
   auto dense_cuda_solver = CUDADenseCholesky::Create(options);
   ASSERT_NE(dense_cuda_solver, nullptr);
@@ -97,6 +103,8 @@ TEST(CUDADenseCholesky, NegativeMatrix) {
         0, 0, -1;
   const Eigen::Vector3d b = Eigen::Vector3d::Ones();
   LinearSolver::Options options;
+  ContextImpl context;
+  options.context = &context;
   options.dense_linear_algebra_library_type = CUDA;
   auto dense_cuda_solver = CUDADenseCholesky::Create(options);
   ASSERT_NE(dense_cuda_solver, nullptr);
@@ -110,6 +118,8 @@ TEST(CUDADenseCholesky, NegativeMatrix) {
 TEST(CUDADenseCholesky, MustFactorizeBeforeSolve) {
   const Eigen::Vector3d b = Eigen::Vector3d::Ones();
   LinearSolver::Options options;
+  ContextImpl context;
+  options.context = &context;
   options.dense_linear_algebra_library_type = CUDA;
   auto dense_cuda_solver = CUDADenseCholesky::Create(options);
   ASSERT_NE(dense_cuda_solver, nullptr);

internal/ceres/cuda_dense_qr_test.cc

@@ -56,6 +56,8 @@ TEST(CUDADenseQR, QR4x4Matrix) {
         0,   0,   0, 1;
   const Eigen::Vector4d b = Eigen::Vector4d::Ones();
   LinearSolver::Options options;
+  ContextImpl context;
+  options.context = &context;
   options.dense_linear_algebra_library_type = CUDA;
   auto dense_cuda_solver = CUDADenseQR::Create(options);
   ASSERT_NE(dense_cuda_solver, nullptr);
@@ -85,6 +87,8 @@ TEST(CUDADenseQR, QR4x2Matrix) {
         0,   0;
   const std::vector<double> b(4, 1.0);
   LinearSolver::Options options;
+  ContextImpl context;
+  options.context = &context;
   options.dense_linear_algebra_library_type = CUDA;
   auto dense_cuda_solver = CUDADenseQR::Create(options);
   ASSERT_NE(dense_cuda_solver, nullptr);
@@ -107,6 +111,8 @@ TEST(CUDADenseQR, QR4x2Matrix) {
 TEST(CUDADenseQR, MustFactorizeBeforeSolve) {
   const Eigen::Vector3d b = Eigen::Vector3d::Ones();
   LinearSolver::Options options;
+  ContextImpl context;
+  options.context = &context;
   options.dense_linear_algebra_library_type = CUDA;
   auto dense_cuda_solver = CUDADenseQR::Create(options);
   ASSERT_NE(dense_cuda_solver, nullptr);
@@ -123,6 +129,8 @@ TEST(CUDADenseQR, Randomized1600x100Tests) {
   using SolutionType = Eigen::Matrix<double, Eigen::Dynamic, 1>;
 
   LinearSolver::Options options;
+  ContextImpl context;
+  options.context = &context;
   options.dense_linear_algebra_library_type = ceres::CUDA;
   std::unique_ptr<DenseQR> dense_qr = CUDADenseQR::Create(options);
 

internal/ceres/dense_cholesky.cc

@@ -36,6 +36,7 @@
 #include <vector>
 
 #ifndef CERES_NO_CUDA
+#include "ceres/context_impl.h"
 #include "cuda_runtime.h"
 #include "cusolverDn.h"
 #endif  // CERES_NO_CUDA
@@ -193,36 +194,17 @@ LinearSolverTerminationType LAPACKDenseCholesky::Solve(const double* rhs,
 
 #ifndef CERES_NO_CUDA
 
-bool CUDADenseCholesky::Init(std::string* message) {
-  if (cusolverDnCreate(&cusolver_handle_) != CUSOLVER_STATUS_SUCCESS) {
-    *message = "cuSolverDN::cusolverDnCreate failed.";
-    return false;
-  }
-  if (cudaStreamCreateWithFlags(&stream_, cudaStreamNonBlocking) !=
-      cudaSuccess) {
-    *message = "CUDA::cudaStreamCreateWithFlags failed.";
-    cusolverDnDestroy(cusolver_handle_);
-    return false;
-  }
-  if (cusolverDnSetStream(cusolver_handle_, stream_) !=
-      CUSOLVER_STATUS_SUCCESS) {
-    *message = "cuSolverDN::cusolverDnSetStream failed.";
-    cudaStreamDestroy(stream_);
-    cusolverDnDestroy(cusolver_handle_);
+bool CUDADenseCholesky::Init(ContextImpl* context, std::string* message) {
+  if (!context->InitCUDA(message)) {
     return false;
   }
+  cusolver_handle_ = context->cusolver_handle_;
+  stream_ = context->stream_;
   error_.Reserve(1);
   *message = "CUDADenseCholesky::Init Success.";
   return true;
 }
 
-CUDADenseCholesky::~CUDADenseCholesky() {
-  if (cusolver_handle_ != nullptr) {
-    CHECK_EQ(cusolverDnDestroy(cusolver_handle_), CUSOLVER_STATUS_SUCCESS);
-    CHECK_EQ(cudaStreamDestroy(stream_), cudaSuccess);
-  }
-}
-
 LinearSolverTerminationType CUDADenseCholesky::Factorize(
     int num_cols, double* lhs, std::string* message) {
   factorize_result_ = LinearSolverTerminationType::LINEAR_SOLVER_FATAL_ERROR;
@@ -326,7 +308,7 @@ std::unique_ptr<CUDADenseCholesky> CUDADenseCholesky::Create(
   auto cuda_dense_cholesky =
       std::unique_ptr<CUDADenseCholesky>(new CUDADenseCholesky());
   std::string cuda_error;
-  if (cuda_dense_cholesky->Init(&cuda_error)) {
+  if (cuda_dense_cholesky->Init(options.context, &cuda_error)) {
     return cuda_dense_cholesky;
   }
   // Initialization failed, destroy the object (done automatically) and return a

internal/ceres/dense_cholesky.h

@@ -44,6 +44,7 @@
 #include "ceres/linear_solver.h"
 #include "glog/logging.h"
 #ifndef CERES_NO_CUDA
+#include "ceres/context_impl.h"
 #include "cuda_runtime.h"
 #include "cusolverDn.h"
 #endif  // CERES_NO_CUDA
@@ -140,7 +141,6 @@ class CERES_NO_EXPORT CUDADenseCholesky : public DenseCholesky {
  public:
   static std::unique_ptr<CUDADenseCholesky> Create(
       const LinearSolver::Options& options);
-  ~CUDADenseCholesky() override;
   CUDADenseCholesky(const CUDADenseCholesky&) = delete;
   CUDADenseCholesky& operator=(const CUDADenseCholesky&) = delete;
   LinearSolverTerminationType Factorize(int num_cols,
@@ -152,11 +152,10 @@ class CERES_NO_EXPORT CUDADenseCholesky : public DenseCholesky {
 
  private:
   CUDADenseCholesky() = default;
-  // Initializes the cuSolverDN context, creates an asynchronous stream, and
-  // associates the stream with cuSolverDN. Returns true iff initialization was
-  // successful, else it returns false and a human-readable error message is
-  // returned.
-  bool Init(std::string* message);
+  // Picks up the cuSolverDN and cuStream handles from the context. If
+  // the context is unable to initialize CUDA, returns false with a
+  // human-readable message indicating the reason.
+  bool Init(ContextImpl* context, std::string* message);
 
   // Handle to the cuSOLVER context.
   cusolverDnHandle_t cusolver_handle_ = nullptr;

internal/ceres/dense_cholesky_test.cc

@@ -65,6 +65,8 @@ TEST_P(DenseCholeskyTest, FactorAndSolve) {
   using VectorType = Eigen::Matrix<Scalar, Eigen::Dynamic, 1>;
 
   LinearSolver::Options options;
+  ContextImpl context;
+  options.context = &context;
   options.dense_linear_algebra_library_type = GetParam();
   std::unique_ptr<DenseCholesky> dense_cholesky =
       DenseCholesky::Create(options);

internal/ceres/dense_qr.cc

@@ -34,6 +34,7 @@
 #include <memory>
 #include <string>
 #ifndef CERES_NO_CUDA
+#include "ceres/context_impl.h"
 #include "cusolverDn.h"
 #include "cublas_v2.h"
 #endif  // CERES_NO_CUDA
@@ -310,53 +311,18 @@ LinearSolverTerminationType LAPACKDenseQR::Solve(const double* rhs,
 
 #ifndef CERES_NO_CUDA
 
-bool CUDADenseQR::Init(std::string* message) {
-  if (cublasCreate(&cublas_handle_) != CUBLAS_STATUS_SUCCESS) {
-    *message = "cuBLAS::cublasCreate failed.";
-    return false;
-  }
-  if (cusolverDnCreate(&cusolver_handle_) != CUSOLVER_STATUS_SUCCESS) {
-    *message = "cuSolverDN::cusolverDnCreate failed.";
-    return false;
-  }
-  if (cudaStreamCreateWithFlags(&stream_, cudaStreamNonBlocking) !=
-      cudaSuccess) {
-    *message = "CUDA::cudaStreamCreateWithFlags failed.";
-    cusolverDnDestroy(cusolver_handle_);
-    cublasDestroy(cublas_handle_);
-    return false;
-  }
-  if (cusolverDnSetStream(cusolver_handle_, stream_) != CUSOLVER_STATUS_SUCCESS) {
-    *message = "cuSolverDN::cusolverDnSetStream failed.";
-    cusolverDnDestroy(cusolver_handle_);
-    cudaStreamDestroy(stream_);
-    cublasDestroy(cublas_handle_);
-    return false;
-  }
-  if (cublasSetStream(cublas_handle_, stream_) != CUBLAS_STATUS_SUCCESS) {
-    *message = "cuBLAS::cublasSetStream failed.";
-    cusolverDnDestroy(cusolver_handle_);
-    cublasDestroy(cublas_handle_);
-    cudaStreamDestroy(stream_);
+bool CUDADenseQR::Init(ContextImpl* context, std::string* message) {
+  if (!context->InitCUDA(message)) {
     return false;
   }
+  cublas_handle_ = context->cublas_handle_;
+  cusolver_handle_ = context->cusolver_handle_;
+  stream_ = context->stream_;
   error_.Reserve(1);
   *message = "CUDADenseQR::Init Success.";
   return true;
 }
 
-CUDADenseQR::~CUDADenseQR() {
-  if (cublas_handle_ != nullptr) {
-    CHECK_EQ(cublasDestroy(cublas_handle_), CUBLAS_STATUS_SUCCESS);
-  }
-  if (cusolver_handle_ != nullptr) {
-    CHECK_EQ(cusolverDnDestroy(cusolver_handle_), CUSOLVER_STATUS_SUCCESS);
-  }
-  if (stream_ != nullptr) {
-    CHECK_EQ(cudaStreamDestroy(stream_), cudaSuccess);
-  }
-}
-
 LinearSolverTerminationType CUDADenseQR::Factorize(
     int num_rows, int num_cols, double* lhs, std::string* message) {
   factorize_result_ = LinearSolverTerminationType::LINEAR_SOLVER_FATAL_ERROR;
@@ -496,7 +462,7 @@ std::unique_ptr<CUDADenseQR> CUDADenseQR::Create(
   auto cuda_dense_qr =
       std::unique_ptr<CUDADenseQR>(new CUDADenseQR());
   std::string cuda_error;
-  if (cuda_dense_qr->Init(&cuda_error)) {
+  if (cuda_dense_qr->Init(options.context, &cuda_error)) {
     return cuda_dense_qr;
   }
   // Initialization failed, destroy the object (done automatically) and return a

internal/ceres/dense_qr.h

@@ -45,7 +45,9 @@
 #include "ceres/internal/export.h"
 #include "ceres/linear_solver.h"
 #include "glog/logging.h"
+
 #ifndef CERES_NO_CUDA
+#include "ceres/context_impl.h"
 #include "ceres/cuda_buffer.h"
 #include "cuda_runtime.h"
 #include "cublas_v2.h"
@@ -153,7 +155,6 @@ class CERES_NO_EXPORT CUDADenseQR : public DenseQR {
  public:
   static std::unique_ptr<CUDADenseQR> Create(
       const LinearSolver::Options& options);
-  ~CUDADenseQR() override;
   CUDADenseQR(const CUDADenseQR&) = delete;
   CUDADenseQR& operator=(const CUDADenseQR&) = delete;
   LinearSolverTerminationType Factorize(int num_rows,
@@ -166,11 +167,10 @@ class CERES_NO_EXPORT CUDADenseQR : public DenseQR {
 
  private:
   CUDADenseQR();
-  // Initializes the cuSolverDN context, creates an asynchronous stream, and
-  // associates the stream with cuSolverDN. Returns true iff initialization was
-  // successful, else it returns false and a human-readable error message is
-  // returned.
-  bool Init(std::string* message);
+  // Picks up the cuSolverDN, cuBLAS, and cuStream handles from the context. If
+  // the context is unable to initialize CUDA, returns false with a
+  // human-readable message indicating the reason.
+  bool Init(ContextImpl* context,std::string* message);
 
   // Handle to the cuSOLVER context.
   cusolverDnHandle_t cusolver_handle_ = nullptr;

internal/ceres/dense_qr_test.cc

@@ -67,6 +67,8 @@ TEST_P(DenseQRTest, FactorAndSolve) {
   using VectorType = Eigen::Matrix<Scalar, Eigen::Dynamic, 1>;
 
   LinearSolver::Options options;
+  ContextImpl context;
+  options.context = &context;
   options.dense_linear_algebra_library_type = GetParam();
   const double kEpsilon = std::numeric_limits<double>::epsilon() * 1.5e4;
   std::unique_ptr<DenseQR> dense_qr = DenseQR::Create(options);