nvidia-cusolver

import cupy as cp from cuda import cusolver, cuda # Initialize CUDA context (often implicit with CuPy) # Create a symmetric positive-definite matrix on GPU using CuPy # For simplicity, let's create a diagonally dominant matrix n = 4 A_host = cp.array([ [4.0, 1.0, 1.0, 1.0], [1.0, 3.0, 0.0, 0.0], [1.0, 0.0, 2.0, 0.0], [1.0, 0.0, 0.0, 1.0] ], dtype=cp.float32) A_device = cp.asarray(A_host, dtype=cp.float32) # Allocate memory for factorization output (in-place for potrf) # Allocate memory for info (error code) info_device = cp.zeros(1, dtype=cp.int32) # Create a cuSOLVER handle handle = None try: handle = cusolver.create_handle() # Query workspace size for potrf (Cholesky factorization) lwork = cusolver.spotrf_bufferSize(handle, cusolver.cudaSolver_fact_info.CUSOLVER_STATUS_SUCCESS, n, A_device.data.ptr, n) workspace = cp.zeros(lwork, dtype=cp.float32) # Perform Cholesky factorization: A = L * L^T (or U^T * U) # We'll use CUSOLVER_FILL_MODE_LOWER (lower triangle) cusolver.spotrf(handle, cusolver.cudaSolver_fact_info.CUSOLVER_FILL_MODE_LOWER, n, A_device.data.ptr, n, workspace.data.ptr, lwork, info_device.data.ptr) # Check for errors info = info_device.get() if info[0] != 0: print(f"Cholesky factorization failed with error code: {info[0]}") else: print("Original Matrix (GPU):\n", A_host) print("Cholesky Factor L (GPU, lower triangle of A_device):\n", A_device.get()) finally: if handle: cusolver.destroy_handle(handle)