CuPy: NumPy and SciPy for GPU(github.com)

    import cupy as cp

    vector_add = cp.RawKernel("""
    extern "C" __global__
    void vector_add(const float *A, const float *B, float *C, int num_elements) {
        int i = blockDim.x * blockIdx.x + threadIdx.x;

        if (i < num_elements) { C[i] = A[i] + B[i]; }
    }
    """, "vector_add")

    num_elements = 50_000

    block_size = 256
    # round up to next multiple of block_size
    grid_size = (num_elements + block_size - 1) // block_size

    a = cp.random.rand(num_elements, dtype=cp.float32)
    b = cp.random.rand(num_elements, dtype=cp.float32)
    c = cp.zeros(num_elements, dtype=cp.float32)

    args = (a, b, c, num_elements)

    print(f"[Vector addition of {num_elements} elements]")
    print(f"CUDA kernel launch with {grid_size} blocks of {block_size} threads each")

    vector_add((grid_size,), (block_size,), args)

    incorrect = cp.abs(a + b - c) > 1e-5

    if cp.any(incorrect):
        print("Result verification failed at element", cp.argmax(incorrect))

    print("Test PASSED")
    print("SUCCESS")

    ```cpp
    <code here>
    ```

    import scipy.linalg as la

    if cuda:
        import cupy as cp
        import cupy.linalg as cla

        ε = cp.asnumpy(cla.eigvalsh(cp.asarray(H)))
    else:
        ε = la.eigvalsh(H)