python tensorflow pytorch pyopengl pyopencl

How to run Python script on a Discrete Graphics AMD GPU?

WHAT I WANT TO DO:

I have a script that I use for factorizing prime numbers given a certain range:

# Python program to display all the prime numbers within an interval

lower = 900
upper = 1000

print("Prime numbers between", lower, "and", upper, "are:")

for num in range(lower, upper + 1):
   # all prime numbers are greater than 1
   if num > 1:
       for i in range(2, num):
           if (num % i) == 0:
               break
       else:
           print(num)

I would like to use the GPU instead of the CPU to run such script so it would be faster

THE PROBLEM:

I don't have a NVIDIA GPU on my Intel NUC NUC8i7HVK but a "Discrete GPU"

If I run this code to check what are my GPUs:

import pyopencl as cl
import numpy as np

a = np.arange(32).astype(np.float32)
res = np.empty_like(a)

ctx = cl.create_some_context()
queue = cl.CommandQueue(ctx)

mf = cl.mem_flags
a_buf = cl.Buffer(ctx, mf.READ_ONLY | mf.COPY_HOST_PTR, hostbuf=a)
dest_buf = cl.Buffer(ctx, mf.WRITE_ONLY, res.nbytes)

prg = cl.Program(ctx, """
    __kernel void sq(__global const float *a,
    __global float *c)
    {
      int gid = get_global_id(0);
      c[gid] = a[gid] * a[gid];
    }
    """).build()

prg.sq(queue, a.shape, None, a_buf, dest_buf)

cl.enqueue_copy(queue, res, dest_buf)

print (a, res)

I receive:

[0] <pyopencl.Platform 'AMD Accelerated Parallel Processing' at 0x7ffb3d492fd0>
[1] <pyopencl.Platform 'Intel(R) OpenCL HD Graphics' at 0x187b648ed80>

THE POSSIBLE APPROACH TO THE PROBLEM:

I found a guide that takes you by the hand and explains step by step how to run it on your GPU. But all Pyhton libraries that pipes Python through the GPU like PyOpenGL, PyOpenCL, Tensorflow (Force python script on GPU), PyTorch, etc... are tailored for NVIDIA.

In case you have an AMD all libraries ask for ROCm but such software still doesn't support integrated GPU or Discrete GPU as far as I know (see my own reply below).

I only found a guide that talks about such approach but I cannot make it work.

Is there hope or I'm just tying to do something impossible?

EDIT: Reply to @chapelo

If I choose 0 the reply is:

Set the environment variable PYOPENCL_CTX='0' to avoid being asked again.
[ 0.  1.  2.  3.  4.  5.  6.  7.  8.  9. 10. 11. 12. 13. 14. 15. 16. 17.
 18. 19. 20. 21. 22. 23. 24. 25. 26. 27. 28. 29. 30. 31.] [  0.   1.   4.   9.  16.  25.  36.  49.  64.  81. 100. 121. 144. 169.
 196. 225. 256. 289. 324. 361. 400. 441. 484. 529. 576. 625. 676. 729.
 784. 841. 900. 961.]

If I choose 1 the reply is:

Set the environment variable PYOPENCL_CTX='1' to avoid being asked again.
[ 0.  1.  2.  3.  4.  5.  6.  7.  8.  9. 10. 11. 12. 13. 14. 15. 16. 17.
 18. 19. 20. 21. 22. 23. 24. 25. 26. 27. 28. 29. 30. 31.] [  0.   1.   4.   9.  16.  25.  36.  49.  64.  81. 100. 121. 144. 169.
 196. 225. 256. 289. 324. 361. 400. 441. 484. 529. 576. 625. 676. 729.
 784. 841. 900. 961.]

Solution

The following code is sample of a complete python program that usually includes:

The import statements
The function definitions
The main() function
The if __name__ == "__main__": section.

I hope this helps you solve your problem.

import pyprimes
from math import sqrt
import numpy as np

import pyopencl as cl
import pyopencl.algorithm
import pyopencl.array

def primes_below(number):
    """Generate a list of prime numbers below a specified  `number`"""
    n = lambda a: 2 if a==0 else 2*a + 1
    limit = int(sqrt(number)) + 1
    size = number//2
    primes = [True] * size
    for i in range(1, size):
        if primes[i]:
            num = n(i)
            if num > limit: break
            for j in range(i+num, size, num):
                primes[j] = False
    for i, flag in enumerate(primes):
        if flag:
            yield n(i)

def primes_between(lo, hi):
    """Generate a list of prime numbers betwenn `lo` and `hi` numbers"""
    primes = list(primes_below(int(sqrt(hi))+1))
    size = (hi - lo - (0 if hi%2 else 1))//2 + 1
    n = lambda a: 2*a + lo + (0 if lo%2 else 1)
    numbers = [True]*size
    for i, prime in enumerate(primes):
        if i == 0: continue # avoid dividing by 2
        nlo = n(0)
        # slower # start = prime * (nlo//prime + 1) if nlo%prime else 0
        start = 0
        while (n(start)%prime) != 0: 
            start += 1
        for j in range(start, size, prime):
            numbers[j] = False
    for i, flag in enumerate(numbers):
        if flag:
            yield n(i)

def primes_between_using_cl(lo, hi):
    """Generate a list of prime numbers betwenn a lo and hi numbers
    this is a parallel algorithm using pyopencl"""
    primes = list(primes_below(int(sqrt(hi))+1))
    size_primes_h = np.array( (len(primes)-1, ), dtype=np.int32)
    numbers_h = np.arange(  lo + (0 if lo&1 else 1), 
                                  hi + (0 if hi&1 else 1),
                                  2,
                                  dtype=np.int32)
    size = (hi - lo - (0 if hi%2 else 1))//2 + 1
    code = """\
    __kernel 
    void is_prime( __global const int *primes,
                        __global         int *numbers) {
      int gid = get_global_id(0);
      int num = numbers[gid];
      int max = (int) (sqrt((float)num) + 1.0);
      for (; *primes; ++primes) {
         if (*primes > max) break;
         if (num % *primes == 0) {
            numbers[gid] = 0;
            return;
         }
      }
    }
    """
    platforms = cl.get_platforms()
    ctx = cl.Context(dev_type=cl.device_type.ALL,
        properties=[(cl.context_properties.PLATFORM, platforms[0])])
    queue = cl.CommandQueue(ctx)
    prg = cl.Program(ctx, code).build()
    numbers_d = cl.array.to_device(queue, numbers_h)
    primes_d = cl.array.to_device(queue, np.array(primes[1:], dtype=np.int32))
    prg.is_prime(queue, (size, ), None, primes_d.data, numbers_d.data)
    array, length = cl.algorithm.copy_if(numbers_d, "ary[i]>0")[:2]
    yield from array.get()[:length.get()]

def test(f, lo, hi):
    """Test that all prime numbers are generated by comparing with the
    output of the library `pyprimes`"""
    a = filter(lambda p: p>lo, pyprimes.primes_below(hi))
    b = f(lo, hi)
    result = True
    for p, q in zip (a, b):
        if p != q:
            print(p, q)
            result = False
    return result
    
def main():
    lower = 1000
    upper = 5000
    print("The prime numbers between {} and {}, are:".format(lower,upper))
    print()
    for p in primes_between_using_cl(lower, upper):
        print(p, end=' ')
    print()

if __name__ == '__main__':
    main()