프랙탈 (Fractal)
프랙탈은 자기 유사성을 가지는 기하학적 구조입니다. 자기 유사성은 같은 패턴에 대해 재귀 또는 반복을 이용해 구현합니다. 아래 그림은 처음 삼각형의 한변에 대해 , 지점에서 새로운 삼각형이 만들어지는 반복적인 패턴을 보여주고 있습니다.
이번 글에서 다룰 프랙탈은 망델브로(Mandelbrot) 집합입니다. 망델브로 집합은 이미지 버퍼의 각각의 픽셀을 복소 평면의 좌표로 변환 후 반복연산을 수행해 그릴 수 있습니다. 망델브로 집합은 픽셀마다 같은 연산을 수행하므로 병렬 프로그래밍을 적용하기에 적합한 예제입니다.
각각의 픽셀에서는 아래와 같이 복소수 과 에 대한 수열 연산이 수행됩니다.
위 식의 상수 의 값으로서 각각의 이미지 버퍼 상 픽셀에 대응하는 복소수를 대입하게 됩니다. 상수 를 점화식의 초기값으로 시작해 의 절댓값이 2보다 크면 수열은 발산한다고 할 수 있습니다. 이 때 2라는 값은 발산하는 수열의 계산을 미리 막아 주는 역할을 하며, 경계값이라고 부릅니다. 수열이 발산하면 는 집합에 속하지 않습니다. 반면에 의 절댓값이 2보다 작을 경우 계산이 무한히 반복되므로 반복 횟수에 제한을 둬야 합니다. 이를 반영한 식은 병렬 연산을 수행할 OpenCL 의 커널에 작성하게 됩니다.
병렬 연산 코드
0. 복소수 구조체 정의
typedef struct s_complex
{
float r;
float i;
} t_complex;구조체 멤버 r 은 실수부, i 는 허수부를 뜻합니다.
1. 복소 평면의 범위 설정
이미지 버퍼가 띄워질 윈도우의 가로 세로 크기를 고려해 복소 평면의 범위를 설정합니다.
# define COORD_WIDTH 5.0f // 복소 평면의 실수부 넓이
# define WINDOW_WIDTH 1000 // 윈도우 넓이
# define WINDOW_HEIGHT 800 // 윈도우 높이
...
t_complex init_entry_point(void)
{
t_complex point;
float coord_height;
float window_width = WINDOW_WIDTH;
point.r = -1.0f * COORD_WIDTH / 2.0f;
coord_height = COORD_WIDTH * WINDOW_HEIGHT / window_width;
point.i = -1.0f * coord_height / 2.0f;
return (point);
}
float init_delta(void)
{
return (COORD_WIDTH / WINDOW_WIDTH);
}함수 init_entry_point 는 매크로 상수로 정의된 윈도우 넓이, 높이와 복소 평면의 실수부 넓이를 바탕으로 복소평면의 가장 왼쪽, 위의 복소수 좌표를 반환합니다. 이는 이미지 버퍼의 가장 왼쪽, 위 픽셀과 대응됩니다. init_delta 는 윈도우 넓이와 복소평면의 넓이의 비율을 이용해 픽셀 사이 간격의 복소 평면 상 비율을 반환합니다. 두 함수의 반환값과 OpenCL 연산 유닛의 global 인덱스를 이용해 각 픽셀에 대응하는 복소수를 구할 수 있습니다.
2. OpenCL 오브젝트 생성
오브젝트 생성은 OpenCL 기초 (1) 과 OpenCL 기초 (2) 에서 연습용 예제를 통해 작성한 코드와 크게 다르지 않습니다. 다만 이미지 버퍼를 윈도우에 출력하는 과정은 이 글에서 생략하겠습니다 (필자의 경우 외부에서 제공하고 있는 그래픽 라이브러리를 사용). 병렬 연산의 결과인 이미지 버퍼가 반환된 후엔 OpenGL 또는 Vulkan 과 같은 익숙한 그래픽 라이브러리를 사용해 윈도우에 출력하시면 됩니다.
이어지는 글에서 다루게 될 OpenCL API 의 자세한 설명은 이전 포스트나 공식문서를 참고하시기 바랍니다.
오브젝트 생성을 위한 첫 단계로 OpenCL API 를 이용해 플랫폼과 디바이스 ID 를 얻습니다. 이번 예제는 GPU 디바이스 하나만 사용합니다.
cl_platform_id platform;
cl_device_id device;
clGetPlatformIDs(1, &platform, NULL);
clGetDeviceIDs(platform, CL_DEVICE_TYPE_GPU, 1, &device, NULL);다음으로 OpenCL 에서 수행되는 커널과 메모리를 관리하는 context 와 실행될 커맨드를 담을 커맨드 큐 구조체를 생성합니다.
int ret;
cl_context context;
cl_command_queue command_queue;
context = clCreateContext(NULL, 1, &device, NULL, NULL, &ret);
command_queue = clCreateCommandQueue(context, device, 0, &ret);그리고 다음 섹션에서 작성하게 될 커널을 사용할 수 있도록 프로그램과 커널 구조체를 생성합니다.
int ret;
cl_program program;
cl_kernel kernel;
char *source_str;
program = clCreateProgramWithSource(context, 1, (const char **)&source_str, NULL, &ret);
clBuildProgram(program, 1, &device, NULL, NULL, NULL);
kernel = clCreateKernel(program, "fractal", &ret);이 때 커널 소스를 문자열 변수 source_str 에 저장하는 과정은 이전 포스트와 마찬가지로 생략하겠습니다. 파일 디스크립터, 파일 포인터를 이용하거나 문자열 상수로서 커널 소스를 작성하셔도 무방합니다. clCreateKernel(program, "fractal", &ret); 에서 "fractal" 은 커널의 이름이고 다음 섹션에서 작성하게 됩니다.
커널의 매개변수이자 병렬 연산의 결과인 메모리 오브젝트를 생성하겠습니다.
int *host_buf = /* 동적 할당 또는 그래픽 라이브러리의 이미지 버퍼 */;
cl_mem dev_mem_obj;
size_t global_work_size = WINDOW_WIDTH * WINDOW_HEIGHT;
dev_mem_obj = clCreateBuffer(context, CL_MEM_WRITE_ONLY,\
sizeof(*host_buf) * global_work_size, NULL, &ret);위 코드로 생성된 메모리 오브젝트는 연산이 끝난 뒤 호스트의 이미지 버퍼 변수 host_buf 에 복사됩니다. global_work_size 는 global 계층의 연산 유닛의 개수를 의미하고 만델브로 집합에 대한 병렬 연산은 각 픽셀에 대해 똑같이 수행되므로 모든 픽셀의 개수인 윈도우의 넓이 높이를 할당합니다.
3. 커널
메모리 오브젝트의 생성까지 병렬 연산의 준비 과정을 마쳤습니다. 이제 만델브로 집합의 복소수 연산 커널을 작성해 보겠습니다.
// fractal.cl
#define MAX_ITERATION 100
#define WHITE 0xFFFFFF
#define BLACK 0x000000
typedef struct s_complex
{
float r;
float i;
} t_complex;커널 소스에 매크로 상수와 구조체를 정의했습니다. OpenCL 커널 코드를 작성하기 위한 언어인 OpenCL C 는 C99 기반의 언어이고 C 에서 지원하는 문법 및 개념을 OpenCL 에서도 사용할 수 있습니다. 구조체 s_complex 는 호스트에서 정의한 s_complex 와 크기, 멤버가 같습니다. 이는 호스트에서 구한 가장 왼쪽 위 픽셀의 복소수 좌표를 데이터의 손실없이 디바이스에 전달하기 위함입니다.
// fractal.cl
...
int iteration(t_complex z, t_complex c)
{
int i;
float aa;
float bb;
float two_ab;
i = 0;
while (i < MAX_ITERATION)
{
aa = z.r * z.r;
bb = z.i * z.i;
two_ab = 2 * z.r * z.i;
if (aa + bb > 4)
return (WHITE);
z.r = aa - bb + c.r;
z.i = two_ab + c.i;
i++;
}
return (BLACK);
}
__kernel void fractal(
__global int *out, // 이미지 버퍼로 쓰일 메모리 오브젝트
t_complex entry_point, // 이미지 버퍼 상 가장 왼쪽 위 픽셀에 대응하는 복소수
float delta, // 픽셀 사이 간격의 복소 평면에 대한 비율
int window_width // 윈도우 넓이
)
{
int idx = get_global_id(0);
int r = idx % window_width;
int i = idx / window_width;
t_complex point;
point.r = entry_point.r + delta * r;
point.i = entry_point.i + delta * i;
out[idx] = iteration(point, point);
}위 코드는 복소수의 만델브로 집합 포함 여부를 판단할 로직을 보여주고 있습니다. 연산 유닛의 개수는 이미지 버퍼의 픽셀 개수와 같습니다. 커널 fractal 에서는 현재 연산 유닛의 global 계층 상 인덱스(int idx = get_global_id(0);) 와 호스트에서 매개 변수로 전달한 윈도우의 넓이(int window_width) 를 이용해 현재 연산 유닛의 복소수, 즉 이미지 버퍼에 대응하는 픽셀의 복소 평면 상 좌표를 구합니다. point.r = ... 과 point.i = ... 코드 구문이 이를 보여주고 있습니다.
함수 iteration 은 커널 fractal 에서 호출하는 함수입니다. OpenCL 에서는 커널이 아닌 이와 같은 함수를 보조함수 (auxiliary function) 라 부릅니다. iteration 은 복소수 연산 를 MAX_ITERATION 번 수행하고 발산하는 수열의 경계값인 2를 이용해 만델브로 집합에 속하지 않는 c 에 대해 흰색 RGB 값을 반환합니다. 조건문 if (aa + bb > 4) 가 이를 보여주고 있습니다. MAX_ITERATION 동안 가 경계값보다 작을 경우 검은색 RGB 값을 반환합니다.
4. 커널 실행
다시 호스트 파트로 돌아와 커널에 필요한 매개변수를 전달할 코드를 작성하겠습니다.
t_complex entry_point = init_entry_point();
float delta = init_delta();
int window_width = WINDOW_WIDTH;
clSetKernelArg(kernel, 0, sizeof(dev_mem_obj), &dev_mem_obj);
clSetKernelArg(kernel, 1, sizeof(entry_point), &entry_point);
clSetKernelArg(kernel, 2, sizeof(delta), &delta);
clSetKernelArg(kernel, 3, sizeof(window_width), &window_width);그리고 커맨드 큐에 커널과 결과 메모리를 읽을 커맨드를 추가합니다.
size_t global_work_size = WINDOW_WIDTH * WINDOW_HEIGHT;
size_t local_work_size = 1;
clEnqueueNDRangeKernel(command_queue, kernel, 1, NULL,\
&global_work_size, &local_work_size, 0, NULL, NULL);
clEnqueueReadBuffer(command_queue, dev_mem_obj, CL_TRUE, 0,\
sizeof(*host_buf) * global_work_size, host_buf, 0, NULL, NULL);커맨드 큐에 추가된 모든 커맨드를 실행합니다.
clFlush(command_queue);
clFinish(command_queue);병렬 연산이 완료된 메모리 오브젝트의 데이터가 이미지 버퍼 host_buf 에 복사되었습니다. 사용하시는 그래픽 라이브러리를 이용해 host_buf 를 윈도우에 출력해보시기 바랍니다. 출력 결과는 아래와 같습니다.
5. 각각의 구조체 해제
마지막으로 각각의 구조체의 메모리를 해제합니다.
clReleaseKernel(kernel);
clReleaseProgram(program);
clReleaseMemObject(dev_mem_obj);
clReleaseCommandQueue(command_queue);
clReleaseContext(context);
아래는 지금까지의 과정을 정리한 코드입니다.
t_complex init_entry_point(void)
{
t_complex point;
float coord_height;
float window_width = WINDOW_WIDTH;
point.r = -1.0f * COORD_WIDTH / 2.0f;
coord_height = COORD_WIDTH * WINDOW_HEIGHT / window_width;
point.i = -1.0f * coord_height / 2.0f;
return (point);
}
float init_delta(void)
{
return (COORD_WIDTH / WINDOW_WIDTH);
}
void example_fractal(void)
{
int fd;
int *host_buf;
char *source_str;
cl_int ret;
cl_platform_id platform;
cl_device_id device;
cl_context context;
cl_command_queue command_queue;
cl_program program;
cl_kernel kernel;
cl_mem dev_mem_obj;
size_t global_work_size = WINDOW_WIDTH * WINDOW_HEIGHT;
size_t local_work_size = 1;
t_complex entry_point;
float delta;
int window_width = WINDOW_WIDTH;
host_buf = /* 동적 할당 또는 그래픽 라이브러리에서 제공하는 버퍼 할당 */;
/*
** example
** host_buf = (int *)malloc(sizeof(int) * global_work_size);
*/
entry_point = init_entry_point();
delta = init_delta();
/* 병렬 연산을 수행할 디바이스 ID 쿼리 */
clGetPlatformIDs(1, &platform, NULL);
clGetDeviceIDs(platform, CL_DEVICE_TYPE_GPU, 1, &device, NULL);
/* context 와 커맨드 큐 구조체 생성 */
context = clCreateContext(NULL, 1, &device, NULL, NULL, &ret);
command_queue = clCreateCommandQueue(context, device, 0, &ret);
/*
** 아래는 .cl 파일의 소스를 문자열로 복사하는 코드입니다.
** get_file_content 와 같은 기능을 하는 함수는 어렵지 않게 만드실 수 있습니다.
*/
fd = open("kernels/fractal.cl", O_RDONLY);
source_str = get_file_content(fd);
close(fd);
/* 프로그램 구조체 생성 */
program = clCreateProgramWithSource(context, 1, (const char **)&source_str, NULL, &ret);
clBuildProgram(program, 1, &device, NULL, NULL, NULL);
/* 커널 구조체 생성 */
kernel = clCreateKernel(program, "fractal", &ret);
/* 메모리 오브젝트 생성 */
dev_mem_obj = clCreateBuffer(context, CL_MEM_WRITE_ONLY,\
sizeof(*host_buf) * global_work_size, NULL, &ret);
/* 커널 매개변수 설정 */
clSetKernelArg(kernel, 0, sizeof(dev_mem_obj), &dev_mem_obj);
clSetKernelArg(kernel, 1, sizeof(entry_point), &entry_point);
clSetKernelArg(kernel, 2, sizeof(delta), &delta);
clSetKernelArg(kernel, 3, sizeof(window_width), &window_width);
/* 커맨드 큐에 커널과 메모리 읽는 커맨드 추가 */
clEnqueueNDRangeKernel(command_queue, kernel, 1, NULL,\
&global_work_size, &local_work_size, 0, NULL, NULL);
clEnqueueReadBuffer(command_queue, dev_mem_obj, CL_TRUE, 0,\
sizeof(*host_buf) * global_work_size, host_buf, 0, NULL, NULL);
/* 커맨드 큐에 추가된 커맨드 실행 */
clFlush(command_queue);
clFinish(command_queue);
/* 윈도우에 host_buf 출력 */
/* 필요한 구조체 메모리 해제 */
clReleaseKernel(kernel);
clReleaseProgram(program);
clReleaseMemObject(dev_mem_obj);
clReleaseCommandQueue(command_queue);
clReleaseContext(context);
}kernels/fractal.cl
#define MAX_ITERATION 100
#define WHITE 0xFFFFFF
#define BLACK 0x000000
typedef struct s_complex
{
float r;
float i;
} t_complex;
int iteration(t_complex z, t_complex c)
{
int i;
float aa;
float bb;
float two_ab;
i = 0;
while (i < MAX_ITERATION)
{
aa = z.r * z.r;
bb = z.i * z.i;
two_ab = 2 * z.r * z.i;
if (aa + bb > 4)
return (WHITE);
z.r = aa - bb + c.r;
z.i = two_ab + c.i;
i++;
}
return (BLACK);
}
__kernel void fractal(
__global int *out, // 이미지 버퍼로 쓰일 메모리 오브젝트
t_complex entry_point, // 이미지 버퍼 상 가장 왼쪽 위 픽셀에 대응하는 복소수
float delta, // 각 픽셀 간의 복소 평면에 대응하는 간격
int window_width // 윈도우 넓이
)
{
int idx = get_global_id(0);
int r = idx % window_width;
int i = idx / window_width;
t_complex point;
point.r = entry_point.r + delta * r;
point.i = entry_point.i + delta * i;
out[idx] = iteration(point, point);
}
