지난 번 글에서는 간단하게 왜 fractal 그리는 작업에 GPU가 필요한지에 대해 간단히 이야기를 해봤습니다. 이번에는 그 fractal을 그리는데 사용되는 CUDA에 대해 설명을 해보도록 하겠습니다.


  CUDA(Compute Unified Device Architecture)는 C언어(정확히는 CUDA C)를 사용해서 GPU에서 알고리즘을 작성할 수 있는 GPGPU 기술입니다. 정확히는 오직 Nvidia GPU에서만 지원하는 기술이죠. AMD GPU나 기타 다른 GPU에서도 비슷한 역할을 하는 OpenCL(Nvidia GPU에서도 돌아갑니다)이라는 것도 있지만, 제가 알기로는 CUDA가 더 좋다고 알고 있습니다. 요즘엔 딥 러닝 때문에 기본적으로 연산을 위해 GPU를 쓰는 경우가 있는데, 대체적으로 딥 러닝 프레임워크가 CUDA를 지원해서 딥 러닝 할 때는 어쩔 수 없이 Nvidia GPU를 고를 수 밖에 없는 일이 발생하기도 했습니다.


  원래 GPU는 CPU와 다르기 때문에 기본적으로 CUDA를 이용해 병렬 프로그래밍을 짤 때에도 기존 프로그래밍과는 조금 다른 방식으로 프로그래밍을 해야 합니다. 무엇보다 RAM을 동시에 같이 쓰지 않기 때문에 CUDA에서는 컴퓨터(Host)와 GPU(Device) 사이의 메모리 복사를 통해 작업물과 결과물을 주고 받습니다. 그리고 Kernel이라는 GPU에서의 function을 짜고 이를 GPU에서 돌려서 연산을 수행합니다. 즉, 기본적인 알고리즘은 아래와 같이 이뤄지게 됩니다.


1. CPU에서 GPU로 데이터 전송(Host to Device 메모리 복사)

2. GPU에서 kernel을 통해 데이터 연산을 수행(kernel 연산)

3. 연산 결과물을 GPU에서 CPU로 데이터 전송(Device to Host 메모리 복사)


이 과정을 원래는 CUDA의 기본 지원 언어인 C/C++에서 짤 때 일일히 구현해야 하는데, 제가 fractal 그리는 프로그램을 짤 때는 PyCUDA를 사용했기 때문에 저 과정을 간단하게 작성할 수 있었습니다. 제 프로그램에서 CUDA 부분만 가져와보면 아래와 같습니다.


kernel 부분

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
32
__global__ void calcMandelbrot(bool *isJulia, double *c, double *matrix, int *result) {
    int idx_x = blockIdx.x * blockDim.x + threadIdx.x;
    int idx_y = blockIdx.y * blockDim.y + threadIdx.y;
    double a = matrix[0+ (float)idx_x * matrix[2];
    double b = matrix[1+ (float)idx_y * matrix[2];
 
    double an = a, bn = b;
    double aan = an * an, bbn = bn * bn;
 
    for (int i = 0; i < 400; i+=5) {
        // iteration 1
        if (isJulia[0]) { // Julia Set
            bn = 2 * an * bn + c[1];
            an = aan - bbn + c[0];
            aan = an * an;
            bbn = bn * bn;
        } else { // Mandelbrot Set
            bn = 2 * an * bn + b;
            an = aan - bbn + a;
            aan = an * an;
            bbn = bn * bn;
        }
 
        if (an * an + bn * bn > 4.0f) {
            result[idx_x * gridDim.y * blockDim.y + idx_y] = i;
            break;
        }
 
        ... 생략 ...
 
    }
}
cs


main 부분

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
def calcMandelbrot(type, start, size, unit, c=0 + 0j):
    """ Calculating and drawing mandelbrot set and julia set
    :param type: fractal type
    :param start: center point in plane
    :param size: image size(width, height)
    :param unit: gap among pixels
    :param c: julia constant
    :return: image 2D array (format: RGBA)
    """
 
    s = time()
 
    # change size to multiple of 16(bigger or equal to prior size). Because of threads per block is 16
    size = ((size + 255/ 256).astype(np.uint32) * 256
 
    result = np.empty((size[0], size[1]), np.int32)
 
    func = cu.get_function("calcMandelbrot")
    func(cuda.In(np.array([type == Fractal.JULIA], np.bool)),
         cuda.In(np.array([c.real, c.imag], np.float64)),
         cuda.In(np.array([start[0- size[0/ 2 * unit,
                           start[1- size[1/ 2 * unit, unit], np.float64)),
         cuda.Out(result),
         block=(16161), grid=(int(size[0/ 16), int(size[1/ 16)))
 
    # Because in image symmetric transformation occurs between x axis and y axis
    result = np.transpose(result)
 
    return array2imgarray(result, cc.m_cyclic_wrwbw_40_90_c42_s25)
cs


main부분을 보시면 16번째 줄에서 결과물을 복사해서 담을 CPU 메모리를 할당하고 있고, 18번째 줄에서 function을 가져와서 19번째 줄에서 호출을 합니다. 원래 C/C++ 같았으면 19번째 줄에서 파라미터 부분(cuda.In, cuda.Out)을 싸그리 할당하고 복사하는 복잡한 과정을 거쳐야 호출할 수 있는데, PyCUDA에서는 간단하게 GPU에 입력만 될 파라미터는 cuda.In으로, GPU에서 출력만 될 파라미터는 cuda.Out으로 처리만 하면 나머지 할당하고 복사하는 과정을 알아서 해줍니다. 만약에 GPU에 입력된 후 다시 출력될 파라미터 같은 경우는 cuda.InOut으로 해주면 간단하게 처리됩니다.


그런데 잘 보시면 func 호출 시 단순히 파라미터만 있는 것이 아니라 block이나 grid와 같은 것이 보입니다. 이것은 GPU에 얼마 만큼의 병렬 처리 단위를 생성해서 계산을 할 것인지를 정하는 것입니다. 제가 겪었던 첫번째 시행착오가 바로 이 부분이었죠.


CUDA에서 병렬 처리 단위는 block과 thread로 나눠집니다. 기본적으로 여러 개의 block이 있고, 이 block 안에 thread가 있는 것입니다. 그리고 각 thread는 연산을 할 때 자신의 위치를 지표로 삼아 데이터를 불러와서 연산을 합니다. 예를 들어 fractal 이미지 같은 경우 픽셀의 (x, y)값마다 thread가 있을테니 각 thread는 자신의 (x, y)에 대한 값을 연산합니다. CUDA에서는 이 (x, y)와 같은 thread index가 3차원까지 되어 있으며 각 차원 별로, 그리고 각 단위(block, thread) 별로 제한값이 있습니다. 이것은 GPU마다 다른데 특히 block 당 thread의 최대 개수는 최근 device 경우에는 1024개로 제한이 되어 있습니다. 이 제한 값은 Linux에서 CUDA 설치 시 제공되는 Sample 중 deviceQuery라는 것을 컴파일해서 실행해보면 알 수 있습니다. 제 경우(GTX 970)엔 아래와 같이 확인할 수 있습니다.


GTX 970 Device Query


이렇게 각 block과 thread 개수를 정해서 GPU에서 kernel을 실행하면 개수에 맞게 thread가 생성이 되어 kernel 함수가 연산을 수행합니다. 이때 자기 thread가 어디 위치에 있는 것인지 알아야 하기 때문에 CUDA에서는 기본적으로 gridDim, blockIdx, blockDim, threadIdx라는 상수를 제공합니다. 제가 이걸 프로그래밍하면서 당연하게 'blockIdx는 block의 위치값, blockDim은 block의 차원일테니, blockDim으로 block의 총 개수를 보면 되겠다!'라는 생각을 했었는데, 이게 아니더라고요. block의 총 개수: gridDim, thread의 총 개수: blockDim 입니다. blockIdx는 block의 index, threadIdx는 thread의 index이고요. 이 내용을 간단하게 그림으로 나타내면 아래와 같습니다. 3D는 그리기 귀찮아서, 2D로 그려봤습니다.


  암튼 이렇게 CUDA의 기본 구조를 파악하고 나서 간단하게 fractal 그리기 프로그램을 만들었습니다. 초기 버전은 단순하게 block 당 thread를 1개로만 해서 만들었습니다. 이때는 뭐 거의 800x800 그리는데 0.3초 정도로 매우 느렸습니다. 그래서 어떻게 해야 하나 싶어서 CUDA sample을 뒤져서 Mandelbrot 예제(옛날에 CUDA sample 뒤지다가 발견했었거든요 ㅇㅇ)를 찾아내서 코드를 비교해봤습니다. 전체적인 구조는 비슷하더라고요. 차이점을 찾아보니 일단 저는 단일 thread만 사용했고, 이상하게 CUDA sample은 loop를 일부 풀어놨습니다(다음 글에서 소개할 기술입니다). 아무래도 thread가 하나인 게 문제인 거 같아서 thread를 여러 개로 늘려봤더니, 속도가 100배 정도(...) 빨라졌습니다. 위에서 말했다시피 thread는 개수가 겨우 1024개밖에 안 되기 때문에 몇 개일 때 최적의 값을 가질지 한 번 검색을 해봤더니, 알고리즘마다 case by case인 것 같더라고요. 제 알고리즘에서 이를 계산해보니 아래와 같은 그래프가 나왔습니다. 네, 정말, thread 한 개 쓰는 건 무식한 짓이라는 것이 너무 자명하게 보이죠. (지금 다시 그래프를 그려보니 100배까진 차이가 안 나네요, 처음 알고리즘을 도대체 어떻게 짰었길래...)

사실상 42 thread 정도 돌리면 이후는 거의 엇비슷한 속도가 나옵니다. 그 이후는 돌릴 때마다 누가 빠른지 조금씩 다르더라고요. 그래서 저는 fractal 그리는 프로그램은 162 thread로 작성했습니다.


이렇게 간단히 fractal 그리는 프로그램을 만들면서 기본적인 최적화를 해봤습니다. 그리고 CUDA sample에 나와있는 의미심장한 loop를 풀어 둔 것을 검색하면서 또 다른 것을 배웠죠. 다음 글에서 다뤄볼 예정입니다.


오랜만에 글을 써 봅니다. 근래에 학업 공부를 하고 그러다 보니까 블로그에 공부한 걸 잘 안 올렸네요.(나름 자료를 찾아보고 틈틈히 보긴 했는데 쩝)


  2018년 겨울방학에 들어서고, 꽤 전부터 생각해오던 CUDA를 이용한 프렉탈 그리기를 직접 프로그래밍했습니다. 중학교 때 프렉탈을 처음 알게 되고 나서, 그 신비함과 아름다움에 빠져서 언젠가는 꼭 프렉탈을 그려내는 프로그래밍을 하겠다고 다짐을 했는데(그 옛날에 julia set의 상수마다 이미지를 만들어서 1TB 하드를 꽉 채우겠다는 생각도 해본 거 같습니다. 뭐 정작 여러 가지 그려보다 보니 이쁜 모양이 나오는 julia set 상수는 그렇게 많진 않았지만요)겨울방학 들어서고 하게 되었네요.


  요즘엔 안드로이드 프로그래밍을 거의 안 하고(앱 유지보수 정도...) 딥 러닝 등을 Python으로 하고 있습니다. 그래서 이번 프렉탈 그리기는 제가 옛날부터 주로 쓰던 Java AWT라던가 그런 걸로 만들지 않았고, Python 3에다가 PyQt5와 PyCUDA를 활용해서 만들어 봤습니다. 아래 이미지들이 제가 지금 제작한(그리고 제작하고 있는) 앱의 스크린샷입니다.

 

(GitHub 링크: https://github.com/AMakeApp/FractalViewer)


  이러한 프렉탈을 그리는데는 사실 CPU를 써도 크게 상관은 없습니다. 단지 GPU를 쓰는 이유는 프렉탈을 그리는데 들어가는 알고리즘이 상당한 연산을 필요로 하기 때문입니다. 이 연산 과정을 알려면 위의 프렉탈들이 어떻게 만들어지는지 알아볼 필요가 있습니다.


  제가 현재 앱에 구현해 둔 프렉탈은 위 두 이미지인 Mandelbrot Set과 Julia Set입니다. 사실 이 둘은 모양이 상이함에도 상당한 상관관계가 있습니다. 왜냐하면 이 둘의 집합을 계산하는 과정이 상당히 유사하기 때문이죠. 그렇기에 좀 더 많이 알려져 있는 Mandelbrot Set에 대해 알아봅시다. Mandelbrot Set은 이름에서 보이듯이 Mandelbrot가 고안한 집합입니다. 사실 많이 알려져 있듯이 저렇게 알록달록 해야지 Mandelbrot Set은 아니고 수학적으로는 아래 이미지처럼 복소평면 상에서 검정색 점들의 집합이 Mandelbrot Set입니다.

Mandelset hires


이 집합은 일반적으로 우리가 집합을 정의할 때와는 조금 다르게 특이하게 정의되어 있는데, 다음 점화식에서 수열이 발산하지 않게 하는 복소수 $c$의 집합으로 정의되어 있습니다.

$$z_0 = 0, z_{n + 1} = z_n^2 + c$$


즉, 이 집합은 무한 집합일테고 복소평면 상의 어느 한 점에 대해 이 점이 Mandelbrot Set에 속하는지 계산하려면 각 복소수에 대해 수많은 반복 연산을 통해 수열의 발산 여부를 확인해야 합니다. 다행히도 $| z_n | > 4$라면 발산한다고 볼 수 있기 때문에 일정한 반복 연산 후 $z_n$의 절대값을 비교하면 좀 더 수월하게 계산 할 수 있습니다.(아마 이 부분이 수학적으로 증명되어 있는 거 같습니다. 이 문서를 참고하세요) 하지만 그렇다고 해도 수많은 반복 연산이 필요한 것은 사실입니다. 초중고등학교에서 배우는 것처럼 단순히 함수꼴로 구할 수 있는 것도 아니고, 각 점에 대해 일일히 노가다를 해야하는 것과 마찬가지니까요. 각 점에 대해서는 제가 대충 해본 결과 200번 정도 이상하면 눈으로 보기에는 프렉탈이 충분히 이쁘게 나옵니다. 커다란 프렉탈 이미지를 생성하려고 2000 x 2000의 프렉탈 이미지를 생성한다고 합시다. 대략 이정도면 200 * 2000 * 2000 = 8억 번 정도의 연산이 필요로 합니다. 일반적으로 한 알고리즘에 대해서 1억 번 연산에 대한 소요시간을 1초로 잡을텐데(요즘 CPU는 더 빨라서 더 많이 연산을 할 수는 있습니다) 프렉탈 이미지를 계산하기 위해서는 8초나 걸린다는 겁니다. 요즘 CPU는 중사양 CPU도 4코어 8쓰레드 정도 가질텐데, 멀티 코어 프로그래밍을 한다고 해도 최소 1초는 걸린다는 이야기가 됩니다.


  이렇게 CPU는 프렉탈 이미지를 계산하는 데는 상당히 오래 걸리기 때문에 GPU의 수많은 코어를 사용해서 더 빠르게 계산을 할 수 있습니다. 비유를 하자면 CPU는 고급 인력 4~8명이 일하는 것이라면 GPU는 초보자 수 천 명이 일을 하는 것과 비슷합니다. 프렉탈 연산이 많은 반복 연산이 필요로 하지만 정작 그 연산 자체가 그렇게 복잡하지 않기 때문에(특히 분기 구문 같은 것이 많지 않기 때문에) 이러한 경우에는 GPU에서 빠르게 계산할 수 있습니다. 이러한 이유로 프렉탈을 빠르게 그려내기 위해서 GPU를 활용했습니다.


  GPU로 그려내기로 결정하고 나서, 제 컴퓨터에 장착되어 있는 GTX970에서 GPGPU 프로그래밍을 하기 위해서 Nvidia에서 지원하는 CUDA를 활용했습니다. 여기서 CUDA는 일반적으로 CUDA C를 사용해서 C/C++과 함께 사용되는데, GUI로 이미지를 보여주면서 편하게 프로그래밍하고 싶었기 때문에 Python 3를 기반으로 CUDA를 사용할 수 있게 만들어 둔 PyCUDA를 사용했습니다. 그리고 유려한 GUI를 보여주는 Qt를 Python에 접목시킨 PyQt5를 사용했고요.


다음 글부터는 이 Fractal을 그리는데 어떻게 CUDA 프로그래밍을 했는지 제 시행착오와 함께 CUDA에 대해 이야기해보겠습니다.

+ Recent posts

티스토리 툴바