Amacım algoritmayı tam anlamıyla anlatmak değil, sadece bir gerçekleme örneği sağlamak. Bu sebeple daha fazla detaylandırmak istemiyorum. Gerçekleme için C++ dili ve görselleştirme için OpenGL kütüphanesini kullandım. Performans arttırımı için dar boğaz olan dinamik bellek ayrımı (dynamic memory allocation) için bir iyileştirme gerçekleştirdim, fakat bu iyileştirmenin çalışabilmesi için maksimum nokta sınırının belirlenmesi gerekiyor, bu bilinmiyorsa daha farklı bir usul geliştirilebilir veya bu iyileştirme kaldırılabilir. Eğer kararlılık (stability) sizin için daha önemli ise bu iyileştirmeyi kaldırmanızı veya gözden geçirmenizi tavsiye ederim, zira bazı durumlarda (nokta yoğunlaşması gibi) hatalarla karşılaştım ayrıca aynı anda birden fazla ağaç üretmesine mani oluyor (no multiple instance).
Not: Nokta tanımında açık kaynak kodlu ve ücretsiz temin edilebilen GTL (https://code.google.com/p/geometry/) kütüphanesinden istifade ettim.
/* * Point.h * * Created on: Oct 24, 2013 * Author: ramazan */ #ifndef POINT_H_ #define POINT_H_ #include <gtl/vec2.hpp> #include <gtl/vec3.hpp> #include "Utility.h" #include <algorithm> using namespace gtl; class Point { public: Point() { color.setValue(0, 0, 1); location.setValue(randomFloat() * SCREEN_WIDTH, randomFloat() * SCREEN_HEIGHT); } void walk() { location += Vec2f(randomSignedFloat() * 30, randomSignedFloat() * 30); location[0] = std::max(location[0], 0.f); location[0] = std::min(location[0], (float)SCREEN_WIDTH); location[1] = std::max(location[1], 0.f); location[1] = std::min(location[1], (float)SCREEN_HEIGHT); } Vec2f location; Vec3f color; }; #endif /* POINT_H_ */
/* * Utility.h * * Created on: Oct 24, 2013 * Author: ramazan */ #ifndef UTILITY_H_ #define UTILITY_H_ #include <stdlib.h> #include <stdio.h> #define SCREEN_WIDTH 800 #define SCREEN_HEIGHT 600 #define THREESHOLD 0.000001 #define PARMEQ(x, y) (fabs(x - y) <= THREESHOLD) static inline float randomFloat() { return (float)rand()/(float)RAND_MAX; } static inline float randomSignedFloat() { const float r = randomFloat(); return r - 0.5f; } #endif /* UTILITY_H_ */
#ifndef QUADTREE_H_ #define QUADTREE_H_ #include <assert.h> #include <vector> #include "Point.h" // You can remove all optimization lines together. #define MAX_POINT_CAPACITY 1000000 // Optimization line! class QuadNode { public: QuadNode() { leaf = true; point = NULL; children[0] = children[1] = children[2] = children[3] = NULL; } ~QuadNode() { for (int i = 0; i < 4; i++) { QuadNode *child = children[i]; if (child != NULL) delete child; } } void insert(const Point &point) { if (leaf && this->point == NULL) { this->point = &point; } else { makeChild(-1, -1); makeChild(-1, 1); makeChild(1, -1); makeChild(1, 1); QuadNode *node = NULL; if (leaf) { const Point *thisPoint = this->point; this->point = NULL; node = findChild(*thisPoint); node->insert(*thisPoint); } node = findChild(point); node->insert(point); leaf = false; } } void query(const Vec2f ¢er, const float radius, std::vector<const Point *> &points) { const float absDX = fabs(center[0] - this->center[0]); const float absDY = fabs(center[1] - this->center[1]); if (absDX <= (halfSize[0] + radius) && absDY <= (halfSize[1] + radius)) { if (point != NULL) { Vec2f disn = point->location - center; const float dist = disn.length(); if (dist <= radius) points.push_back(point); } else { for (int i = 0; i < 4; i++) { QuadNode *child = children[i]; if (child != NULL) child->query(center, radius, points); } } } } void render() { const Vec2f p0 = center + Vec2f(-halfSize[0], halfSize[1]); const Vec2f p1 = center + Vec2f(-halfSize[0], -halfSize[1]); const Vec2f p2 = center + Vec2f(halfSize[0], -halfSize[1]); const Vec2f p3 = center + Vec2f(halfSize[0], halfSize[1]); glColor3f(1, 0, 0); glBegin(GL_LINES); glVertex2fv(p0.getValue()); glVertex2fv(p1.getValue()); glVertex2fv(p1.getValue()); glVertex2fv(p2.getValue()); glVertex2fv(p2.getValue()); glVertex2fv(p3.getValue()); glVertex2fv(p3.getValue()); glVertex2fv(p0.getValue()); glEnd(); for (int i = 0; i < 4; i++) { QuadNode *child = children[i]; if (child != NULL) child->render(); } } private: int makeChild(const float dX, const float dY) { int index = -1; float newX = halfSize[0] / 2; float newY = halfSize[1] / 2; if (dX < 0) { newX = -newX; if (dY < 0) { index = 1; newY = -newY; } else index = 0; } else { if (dY < 0) { index = 2; newY = -newY; } else index = 3; } verifyChild(index, center + Vec2f(newX, newY)); return index; } QuadNode *findChild(const Point &point) { const float dX = point.location[0] - center[0]; const float dY = point.location[1] - center[1]; return children[makeChild(dX, dY)]; } void verifyChild(const int index, const Vec2f newCenter) { assert(index >= 0 && index <= 3); if (children[index] == NULL) { children[index] = new QuadNode(); assert(children[index] != NULL); children[index]->center = newCenter; children[index]->halfSize = halfSize / 2; } } private: static char memory[]; // Optimization line! static int memoryIndex; // Optimization line! public: void *operator new(size_t size) { // Optimization line! assert(memoryIndex < MAX_POINT_CAPACITY); // Optimization line! return &memory[memoryIndex++ * size]; // Optimization line! } // Optimization line! void operator delete(void *p) { // Optimization line! p = NULL; // Optimization line! } // Optimization line! public: // 0 3 // 1 2 QuadNode *children[4]; const Point *point; Vec2f center; Vec2f halfSize; bool leaf; }; int QuadNode::memoryIndex = 0; // Optimization line! char QuadNode::memory[sizeof(QuadNode) * MAX_POINT_CAPACITY] = { }; // Optimization line! class Quadtree { public: Quadtree(const Vec2f ¢er, const float w, const float h) { root = new QuadNode(); root->center = center; root->halfSize.setValue(w / 2, h / 2); } ~Quadtree() { if (root != NULL) { delete root; } } void insert(const Point &point) { assert(root != NULL); root->insert(point); } void render() { assert(root != NULL); root->render(); } void query(const Vec2f ¢er, const float radius, std::vector<const Point *> &points) { assert(root != NULL); root->query(center, radius, points); } QuadNode *root; }; #endif /* QUADTREE_H_ */
// main.cpp #include <iostream> #include <glut\\glut.h> #include <time.h> #include "Quadtree.h" #define RENDER_WAIT_TIME 25 #define NUMBER_OF_POINTS 30000 using namespace std; Point *points = NULL; Quadtree *qt; Vec2f qCenter(300, 400); float radius = 100; void display(void) { glClear(GL_COLOR_BUFFER_BIT); glLoadIdentity(); glBegin(GL_POINTS); for (int i = 0; i < NUMBER_OF_POINTS; i++) { glColor3fv(points[i].color.getValue()); glVertex2fv(points[i].location.getValue()); } glEnd(); glColor3f(0, 1, 0); glBegin(GL_LINES); for (int i = 0; i < 360; i++) { Vec2f p1 = qCenter + Vec2f(cos((i * M_PI / 180.f)) * radius, sin((i * M_PI / 180.f)) * radius); Vec2f p2 = qCenter + Vec2f(cos(((i + 1) * M_PI / 180.f)) * radius, sin(((i + 1) * M_PI / 180.f)) * radius); glVertex2fv(p1.getValue()); glVertex2fv(p2.getValue()); } glEnd(); if (qt != NULL) qt->render(); //glFlush(); glutSwapBuffers(); } void mouse(int btn, int state, int x, int y) { //if(btn==GLUT_LEFT_BUTTON && state == GLUT_DOWN) axis = 0; } void search() { if (qt == NULL) return; std::vector<const Point *> qpoints; qt->query(qCenter, radius, qpoints); std::cout<<qpoints.size()<<std::endl; for (int i = 0; i < NUMBER_OF_POINTS; i++) { points[i].color.setValue(0, 0, 1); for (unsigned int j = 0; j < qpoints.size(); j++) { if (&points[i] == qpoints[j]) points[i].color.setValue(0, 1, 0); } } } void mainLoop(int value) { glutPostRedisplay(); glutTimerFunc(RENDER_WAIT_TIME, mainLoop, 0); } void build() { qt = new Quadtree(Vec2f(SCREEN_WIDTH / 2, SCREEN_HEIGHT / 2), SCREEN_WIDTH, SCREEN_HEIGHT); for (int i = 0; i < NUMBER_OF_POINTS; i++) qt->insert(points[i]); } void keyPressed(unsigned char key, int x, int y) { if (key == 'x') exit(0); if (key == 'w') qCenter += Vec2f(0, 10); if (key == 's') qCenter += Vec2f(0, -10); if (key == 'a') qCenter += Vec2f(-10, 0); if (key == 'd') qCenter += Vec2f(10, 0); if (key == 'r') { for (int i = 0; i < NUMBER_OF_POINTS; i++) points[i].walk(); build(); } search(); } void myReshape(int w, int h) { glViewport(0, 0, w, h); glMatrixMode(GL_PROJECTION); glLoadIdentity(); glOrtho(0.0f, 800, 0, 600, 1.0f, -1.0f); glMatrixMode(GL_MODELVIEW); } void keySpecial(int key, int x, int y) { if (key == GLUT_KEY_F4) exit(0); } void init() { glShadeModel(GL_FLAT); glPixelStorei(GL_UNPACK_ALIGNMENT, 1); glPointSize(3); points = new Point[NUMBER_OF_POINTS]; //points[0].location.setValue(100, 100); //points[1].location.setValue(110, 120); //points[2].location.setValue(105, 105); build(); search(); } int main(int argc, char** argv) { srand(time(NULL)); glutInit(&argc, argv); glutInitDisplayMode(GLUT_DOUBLE | GLUT_RGBA); glutInitWindowSize(SCREEN_WIDTH, SCREEN_HEIGHT); glutCreateWindow("Quadtree"); glutReshapeFunc(myReshape); glutDisplayFunc(display); glutKeyboardFunc(keyPressed); glutSpecialFunc(keySpecial); init(); glutTimerFunc(RENDER_WAIT_TIME, mainLoop, 0); glutMouseFunc(mouse); glutMainLoop(); return 0; }Bir kaç deneme yapalım:
Bu bölümde 2 boyutlu düzlem üzerindeki geometrik problemlerin çözümünde kullanılan arama işlemini hızlandıran bir veriyapısını inceledik. Aynı veri yapısı 3 boyutlu uzayda da uygulanabilir (Octree - http://en.wikipedia.org/wiki/Octree). Mantık olarak aynı yalnızca 4 yerine 8 eşit parçaya ayrılmaktadır. Performans değerlendirmesini zaman cinsinden değerleri belirterek yapmak istemiyorum, kullandığım bilgisayar üzerinde brute-force yöntemine göre tahmin edebileceğiniz gibi (beklendiği üzere) çok büyük oranda hızlanma sağladı. Bunu kendi probleminize uygulayarak görebilirsiniz. Basit düzeyde anlatmaya çalıştım, umarım işinize yarar.