diff --git a/README.md b/README.md index ee39093..dab7ca4 100644 --- a/README.md +++ b/README.md @@ -1,11 +1,71 @@ -**University of Pennsylvania, CIS 5650: GPU Programming and Architecture, -Project 1 - Flocking** +# University of Pennsylvania, CIS 5650: GPU Programming and Architecture +## Project 1 - Flocking -* (TODO) YOUR NAME HERE - * (TODO) [LinkedIn](), [personal website](), [twitter](), etc. -* Tested on: (TODO) Windows 22, i7-2222 @ 2.22GHz 22GB, GTX 222 222MB (Moore 2222 Lab) +* Zwe Tun + * LinkedIn: https://www.linkedin.com/in/zwe-tun-6b7191256/ +* Tested on: Intel(R) i7-14700HX, 2100 Mhz, RTX 5060 Laptop +![CUDA Flocking](images/CUDA-Flocking.gif) -### (TODO: Your README) +*100,000 boids with Coherent Grid Search* +## Overview +Boids are artificial agents that simulate the behavior of flocking animals. Introduced by Craig Reynolds in 1986, each boid follows three simple rules: +Cohesion - boids move towards the perceived center of mass of their neighbors +Separation - boids avoid getting to close to their neighbors +Alignment - boids generally try to move with the same direction and speed as their neighbors. -Include screenshots, analysis, etc. (Remember, this is public, so don't put -anything here that you don't want to share with the world.) +## Implementation + +### Naive +The naive implementation uses a straightforward algorithm: each boid iterates over every other boid in the system to calculate its updated velocity and position based on the three flocking rules (separation, alignment, and cohesion). +While conceptually simple, this approach results in O(n²) time complexity, making it highly inefficient for large numbers of boids. Every boid must check all others regardless of distance, leading to significant computational overhead. +![CUDA Flocking](images/Naive-CUDA-Flocking.gif) + +*10,000 boids with Naive* + +### Uniform Grid Search +A more effcient algorithm is dividing into 3D cells, and each boid is assigned to a cell based on its position. By storing cell indices and mapping boid data accordingly, each thread can now limit its neighbor search to only nearby cells rather than the entire boid population. +In the scattered version, boid data (position, velocity) is stored in separate buffers, and additional lookup is required to gather information based on the grid cell. This greatly reduces the number of comparisons per boid, improving performance. +![CUDA Flocking](images/Uniform-CUDA-Flocking.gif) + +*10,000 boids with Uniform Grid Search* + +### Coherent Grid Search +The coherent grid improves on the uniform grid approach by taking advantage of spatial locality to optimize memory access on the GPU. In the scattered grid version, boid data is stored in separate buffers and accessed via indirect lookups, causing threads to read from scattered memory locations. This leads to reduced cache utilization. The coherent grid reorganizes the boid data so that boids located in the same or neighboring grid cells are stored contiguously in memory. By reshuffling the boid data arrays to align with their cell indices, the GPU can access data in adjacent memory locations, further improving performance. +![CUDA Flocking](images/Coherant-CUDA-Flocking.gif) + +*10,000 boids with Coherent Grid Search* + +## Performance Analysis + +All performance tests were conducted on Windows 11, Intel(R) i7-14700HX CPU (2.1 GHz), NVIDIA RTX 5060 Laptop GPU. The primary metric used to evaluate performance is Frames Per Second (FPS). Higher FPS values indicate better performance and smoother real-time simulation. + +![CUDA Flocking](images/Off.png) + +![CUDA Flocking](images/On.png) + +![CUDA Flocking](images/Blocks.png) + +![CUDA Flocking](images/CellSearched.png) + +## Questions +### For each implementation, how does changing the number of boids affect performance? Why do you think this is? + +As the number of boids increases, all 3 methods see a decrease in performance (FPS) due to more computational load. This is shown in the plots where all 3 methods see a decrease in FPS as boid counts increases. + +### For each implementation, how does changing the block count and block size affect performance? Why do you think this is? + +Yes, increasing or decreasing the block size tends to slightly lower performance across the board. This happens because the GPU scheduler tries to optimize its task based on the given resources, and if the sizes or counts or not well balanced it can lead to poor optimization. + +### For the coherent uniform grid: did you experience any performance improvements with the more coherent uniform grid? Was this the outcome you expected? Why or why not? + +For the coherent uniform grid, clear performance improvements are found in the graphs above. With 5,000 boids, the coherent grid achieves nearly 2.5× the FPS of the naive implementation. This performance gap is more noticible at higher boid counts—for example, at 100,000 boids, the coherent grid runs at about 55× the FPS of the naive method. This superior performance is expected, as the coherent grid uses spatial locality to get more efficient memory access and reduce unnecessary computations. + +### Did changing cell width and checking 27 vs 8 neighboring cells affect performance? Why or why not? + +I found that checking 8 neighboring cells results in about a 2× increase in FPS compared to checking 27 cells, for both the uniform and coherent grid implementations. The naive implementation’s performance remains unchanged. This improvement is contributed by checking fewer neighboring cells reducing the number of boids each thread must evaluate. In addition checking a smaller area means that neighbor boids that will not affect the overall are skipped wasting less performance. + +## For the Curious: 1,000,000 Boids + +![CUDA Flocking](images/Million-CUDA-Flocking.gif) + +*1,000,000 boids with Coherent Grid Search* diff --git a/images/Block.png b/images/Block.png new file mode 100644 index 0000000..6df08a6 Binary files /dev/null and b/images/Block.png differ diff --git a/images/Blocks.png b/images/Blocks.png new file mode 100644 index 0000000..dc09b35 Binary files /dev/null and b/images/Blocks.png differ diff --git a/images/CUDA-Flocking.gif b/images/CUDA-Flocking.gif new file mode 100644 index 0000000..994ba70 Binary files /dev/null and b/images/CUDA-Flocking.gif differ diff --git a/images/CellSearched.png b/images/CellSearched.png new file mode 100644 index 0000000..2b3df25 Binary files /dev/null and b/images/CellSearched.png differ diff --git a/images/Coherant-CUDA-Flocking.gif b/images/Coherant-CUDA-Flocking.gif new file mode 100644 index 0000000..95641c8 Binary files /dev/null and b/images/Coherant-CUDA-Flocking.gif differ diff --git a/images/Million-CUDA-Flocking.gif b/images/Million-CUDA-Flocking.gif new file mode 100644 index 0000000..bcdb8e5 Binary files /dev/null and b/images/Million-CUDA-Flocking.gif differ diff --git a/images/Naive-CUDA-Flocking.gif b/images/Naive-CUDA-Flocking.gif new file mode 100644 index 0000000..e447bc1 Binary files /dev/null and b/images/Naive-CUDA-Flocking.gif differ diff --git a/images/Off.png b/images/Off.png new file mode 100644 index 0000000..58e7f05 Binary files /dev/null and b/images/Off.png differ diff --git a/images/On.png b/images/On.png new file mode 100644 index 0000000..e7e5a51 Binary files /dev/null and b/images/On.png differ diff --git a/images/Performance.png b/images/Performance.png new file mode 100644 index 0000000..94f2abf Binary files /dev/null and b/images/Performance.png differ diff --git a/images/Uniform-CUDA-Flocking.gif b/images/Uniform-CUDA-Flocking.gif new file mode 100644 index 0000000..642849a Binary files /dev/null and b/images/Uniform-CUDA-Flocking.gif differ diff --git a/src/kernel.cu b/src/kernel.cu index 7149917..0a70f4c 100644 --- a/src/kernel.cu +++ b/src/kernel.cu @@ -14,6 +14,8 @@ #include #include +#include + #include // LOOK-2.1 potentially useful for doing grid-based neighbor search @@ -51,11 +53,11 @@ void checkCUDAError(const char *msg, int line = -1) { // LOOK-1.2 Parameters for the boids algorithm. // These worked well in our reference implementation. -#define rule1Distance 5.0f +#define rule1Distance 7.2f #define rule2Distance 3.0f #define rule3Distance 5.0f -#define rule1Scale 0.01f +#define rule1Scale 0.04f #define rule2Scale 0.1f #define rule3Scale 0.1f @@ -67,10 +69,10 @@ void checkCUDAError(const char *msg, int line = -1) { /*********************************************** * Kernel state (pointers are device pointers) * ***********************************************/ - int numObjects; dim3 threadsPerBlock(blockSize); + // LOOK-1.2 - These buffers are here to hold all your boid information. // These get allocated for you in Boids::initSimulation. // Consider why you would need two velocity buffers in a simulation where each @@ -93,8 +95,11 @@ thrust::device_ptr dev_thrust_particleGridIndices; int *dev_gridCellStartIndices; // What part of dev_particleArrayIndices belongs int *dev_gridCellEndIndices; // to this cell? -// TODO-2.3 - consider what additional buffers you might need to reshuffle +// additional buffers needed to reshuffle // the position and velocity data to be coherent within cells. +glm::vec3 *dev_reshuffledPos; +glm::vec3 *dev_reshuffledVel1; + // LOOK-2.1 - Grid parameters based on simulation parameters. // These are automatically computed for you in Boids::initSimulation @@ -143,11 +148,13 @@ __global__ void kernGenerateRandomPosArray(int time, int N, glm::vec3 * arr, flo } } + /** * Initialize memory, update some globals */ void Boids::initSimulation(int N) { numObjects = N; + dim3 fullBlocksPerGrid((N + blockSize - 1) / blockSize); // LOOK-1.2 - This is basic CUDA memory management and error checking. @@ -166,6 +173,7 @@ void Boids::initSimulation(int N) { dev_pos, scene_scale); checkCUDAErrorWithLine("kernGenerateRandomPosArray failed!"); + // LOOK-2.1 computing grid params gridCellWidth = 2.0f * std::max(std::max(rule1Distance, rule2Distance), rule3Distance); int halfSideCount = (int)(scene_scale / gridCellWidth) + 1; @@ -178,7 +186,29 @@ void Boids::initSimulation(int N) { gridMinimum.y -= halfGridWidth; gridMinimum.z -= halfGridWidth; - // TODO-2.1 TODO-2.3 - Allocate additional buffers here. + // Allocate additional buffers. + cudaMalloc((void**)&dev_particleArrayIndices, N * sizeof(int)); + checkCUDAErrorWithLine("cudaMalloc dev_particleArrayIndices failed!"); + + cudaMalloc((void**)&dev_particleGridIndices, N * sizeof(int)); + checkCUDAErrorWithLine("cudaMalloc dev_particleGridIndices failed!"); + + cudaMalloc((void**)&dev_gridCellStartIndices, gridCellCount * sizeof(int)); + checkCUDAErrorWithLine("cudaMalloc dev_gridCellStartIndices failed!"); + + cudaMalloc((void**)&dev_gridCellEndIndices, gridCellCount * sizeof(int)); + checkCUDAErrorWithLine("cudaMalloc dev_gridCellEndIndices failed!"); + + dev_thrust_particleGridIndices = thrust::device_pointer_cast(dev_particleGridIndices); + dev_thrust_particleArrayIndices = thrust::device_pointer_cast(dev_particleArrayIndices);; + + cudaMalloc((void**)&dev_reshuffledPos, N * sizeof(glm::vec3)); + checkCUDAErrorWithLine("cudaMalloc dev_reshuffledPos failed!"); + + cudaMalloc((void**)&dev_reshuffledVel1 , N * sizeof(glm::vec3)); + checkCUDAErrorWithLine("cudaMalloc dev_vel1 failed!"); + + cudaDeviceSynchronize(); } @@ -243,11 +273,58 @@ __device__ glm::vec3 computeVelocityChange(int N, int iSelf, const glm::vec3 *po // Rule 1: boids fly towards their local perceived center of mass, which excludes themselves // Rule 2: boids try to stay a distance d away from each other // Rule 3: boids try to match the speed of surrounding boids - return glm::vec3(0.0f, 0.0f, 0.0f); + + glm::vec3 velChange(0.0f, 0.0f, 0.0f); + glm::vec3 averagePosition(0.0f, 0.0f, 0.0f); + glm::vec3 averageSpeed(0.0f, 0.0f, 0.0f); + glm::vec3 apart(0.0f, 0.0f, 0.0f); + int rule1neighbors = 0; + int rule3neighbors = 0; + + for (int i = 0; i < N; i++) { + if (i == iSelf) { + continue; + } + + float distance = glm::distance(pos[i], pos[iSelf]); + + if (distance < rule1Distance) { + averagePosition += pos[i]; + rule1neighbors++; + } + + if (distance < rule2Distance) { + apart -= (pos[i] - pos[iSelf]); + } + + if (distance < rule3Distance) { + averageSpeed += vel[i]; + rule3neighbors++; + } + + + } + if (rule1neighbors > 0) { + averagePosition = (averagePosition) / static_cast(rule1neighbors); + velChange += (averagePosition - pos[iSelf]) * rule1Scale; + } + + + velChange += (apart) * rule2Scale; + + + + if (rule3neighbors > 0) { + averageSpeed = (averageSpeed) / static_cast(rule3neighbors); + velChange += (averageSpeed) * rule3Scale; + } + + return velChange; + } /** -* TODO-1.2 implement basic flocking +* basic flocking * For each of the `N` bodies, update its position based on its current velocity. */ __global__ void kernUpdateVelocityBruteForce(int N, glm::vec3 *pos, @@ -255,6 +332,19 @@ __global__ void kernUpdateVelocityBruteForce(int N, glm::vec3 *pos, // Compute a new velocity based on pos and vel1 // Clamp the speed // Record the new velocity into vel2. Question: why NOT vel1? + glm::vec3 velChange; + int index = threadIdx.x + (blockIdx.x * blockDim.x); + if (index >= N) { + return; + } + velChange = computeVelocityChange(N, index, pos, vel1); + glm::vec3 newVel = vel1[index] + velChange; + newVel = glm::clamp(newVel, maxSpeed * -1.0f, maxSpeed); + vel2[index] = newVel; + + + + } /** @@ -295,10 +385,24 @@ __device__ int gridIndex3Dto1D(int x, int y, int z, int gridResolution) { __global__ void kernComputeIndices(int N, int gridResolution, glm::vec3 gridMin, float inverseCellWidth, glm::vec3 *pos, int *indices, int *gridIndices) { - // TODO-2.1 // - Label each boid with the index of its grid cell. // - Set up a parallel array of integer indices as pointers to the actual // boid data in pos and vel1/vel2 + + int index = (blockIdx.x * blockDim.x) + threadIdx.x; + if (index < N) { + int iX = (int)floor((pos[index].x - gridMin.x) * inverseCellWidth); + int iY = (int)floor((pos[index].y - gridMin.y) * inverseCellWidth); + int iZ = (int)floor((pos[index].z - gridMin.z) * inverseCellWidth); + iX = imin(imax(iX, 0), gridResolution - 1); + iY = imin(imax(iY, 0), gridResolution - 1); + iZ = imin(imax(iZ, 0), gridResolution - 1); + int gridIndex = gridIndex3Dto1D(iX, iY, iZ, gridResolution); + indices[index] = index; + gridIndices[index] = gridIndex; + } + + } // LOOK-2.1 Consider how this could be useful for indicating that a cell @@ -312,10 +416,24 @@ __global__ void kernResetIntBuffer(int N, int *intBuffer, int value) { __global__ void kernIdentifyCellStartEnd(int N, int *particleGridIndices, int *gridCellStartIndices, int *gridCellEndIndices) { - // TODO-2.1 + // Identify the start point of each cell in the gridIndices array. // This is basically a parallel unrolling of a loop that goes // "this index doesn't match the one before it, must be a new cell!" + + //Cell start inclusive, cell end exclusive + int index = (blockIdx.x * blockDim.x) + threadIdx.x; + if (index < N) { + if (index == 0) { + gridCellStartIndices[particleGridIndices[index]] = index; + } else if (particleGridIndices[index] != particleGridIndices[index - 1]) { + gridCellStartIndices[particleGridIndices[index]] = index; + gridCellEndIndices[particleGridIndices[index - 1]] = index; + } else if (index == N - 1) { + gridCellEndIndices[particleGridIndices[index]] = index + 1; + } + + } } __global__ void kernUpdateVelNeighborSearchScattered( @@ -324,7 +442,7 @@ __global__ void kernUpdateVelNeighborSearchScattered( int *gridCellStartIndices, int *gridCellEndIndices, int *particleArrayIndices, glm::vec3 *pos, glm::vec3 *vel1, glm::vec3 *vel2) { - // TODO-2.1 - Update a boid's velocity using the uniform grid to reduce + // Update a boid's velocity using the uniform grid to reduce // the number of boids that need to be checked. // - Identify the grid cell that this particle is in // - Identify which cells may contain neighbors. This isn't always 8. @@ -332,17 +450,121 @@ __global__ void kernUpdateVelNeighborSearchScattered( // - Access each boid in the cell and compute velocity change from // the boids rules, if this boid is within the neighborhood distance. // - Clamp the speed change before putting the new speed in vel2 + int index = (blockIdx.x * blockDim.x) + threadIdx.x; + if (index >= N) { + return; + } + + + int iX = (int)floor((pos[index].x - gridMin.x) * inverseCellWidth); + int iY = (int)floor((pos[index].y - gridMin.y) * inverseCellWidth); + int iZ = (int)floor((pos[index].z - gridMin.z) * inverseCellWidth); + iX = imin(imax(iX, 0), gridResolution - 1); + iY = imin(imax(iY, 0), gridResolution - 1); + iZ = imin(imax(iZ, 0), gridResolution - 1); + + glm::vec3 velChange(0.0f, 0.0f, 0.0f); + glm::vec3 averagePosition(0.0f, 0.0f, 0.0f); + glm::vec3 averageSpeed(0.0f, 0.0f, 0.0f); + glm::vec3 apart(0.0f, 0.0f, 0.0f); + int rule1neighbors = 0; + int rule3neighbors = 0; + + + //8 possible cells to check + //loop over z, y + for (int z = -1; z <= 1; z++) { + for (int y = -1; y <= 1; y++) { + + int neighborX = iX; + int neighborY = iY + y; + int neighborZ = iZ + z; + if (neighborX < 0 || neighborY < 0 || neighborZ < 0 || neighborX >= gridResolution || neighborY >= gridResolution || neighborZ >= gridResolution) { + continue; + } + + int gridIndex = gridIndex3Dto1D(neighborX, neighborY, neighborZ, gridResolution); + int startBoid = gridCellStartIndices[gridIndex]; + int endBoid = gridCellEndIndices[gridIndex]; + if (startBoid == -1 || endBoid == -1) { + continue; + } + + + for (int i = startBoid; i < endBoid; i++) { + int neighborBoidIndex = particleArrayIndices[i]; + glm::vec3 neighborPos = pos[neighborBoidIndex]; + glm::vec3 neighborVel1 = vel1[neighborBoidIndex]; + + if (neighborBoidIndex == index) { + continue; + } + + float distance = glm::distance(neighborPos, pos[index]); + + if (distance < rule1Distance) { + averagePosition += neighborPos; + rule1neighbors++; + } + + if (distance < rule2Distance) { + apart -= (neighborPos - pos[index]); + } + + if (distance < rule3Distance) { + averageSpeed += neighborVel1; + rule3neighbors++; + } + + + + } + + + } + } + + + if (rule1neighbors > 0) { + averagePosition /= static_cast(rule1neighbors); + velChange += (averagePosition - pos[index]) * rule1Scale; + } + + + velChange += (apart*rule2Scale); + + if (rule3neighbors > 0) { + averageSpeed /= static_cast(rule3neighbors); + velChange += (averageSpeed*rule3Scale); + } + + // clamp the speed change and record the new velocity into vel2 + glm::vec3 newVel = vel1[index] + velChange; + newVel = glm::clamp(newVel, maxSpeed * -1.0f, maxSpeed); + vel2[index] = newVel; + + } +//reshuffle the position and velocity data to be coherent within cells. +//Memory will now be continuous improving memory access pattern +__global__ void kernReshuffleBoidData(int N, int* particleArrayIndices, + glm::vec3* pos, glm::vec3* vel1, + glm::vec3* reshuffledPos, glm::vec3* reshuffledVel1) { + int index = (blockIdx.x * blockDim.x) + threadIdx.x; + if (index < N) { + int boidIndex = particleArrayIndices[index]; + reshuffledPos[index] = pos[boidIndex]; + reshuffledVel1[index] = vel1[boidIndex]; + } +} + + __global__ void kernUpdateVelNeighborSearchCoherent( int N, int gridResolution, glm::vec3 gridMin, float inverseCellWidth, float cellWidth, int *gridCellStartIndices, int *gridCellEndIndices, glm::vec3 *pos, glm::vec3 *vel1, glm::vec3 *vel2) { - // TODO-2.3 - This should be very similar to kernUpdateVelNeighborSearchScattered, - // except with one less level of indirection. - // This should expect gridCellStartIndices and gridCellEndIndices to refer - // directly to pos and vel1. // - Identify the grid cell that this particle is in // - Identify which cells may contain neighbors. This isn't always 8. // - For each cell, read the start/end indices in the boid pointer array. @@ -351,18 +573,125 @@ __global__ void kernUpdateVelNeighborSearchCoherent( // - Access each boid in the cell and compute velocity change from // the boids rules, if this boid is within the neighborhood distance. // - Clamp the speed change before putting the new speed in vel2 + + int index = (blockIdx.x * blockDim.x) + threadIdx.x; + if (index >= N) { + return; + } + + + int iX = (int)floor((pos[index].x - gridMin.x) * inverseCellWidth); + int iY = (int)floor((pos[index].y - gridMin.y) * inverseCellWidth); + int iZ = (int)floor((pos[index].z - gridMin.z) * inverseCellWidth); + iX = imin(imax(iX, 0), gridResolution - 1); + iY = imin(imax(iY, 0), gridResolution - 1); + iZ = imin(imax(iZ, 0), gridResolution - 1); + + glm::vec3 velChange(0.0f, 0.0f, 0.0f); + glm::vec3 averagePosition(0.0f, 0.0f, 0.0f); + glm::vec3 averageSpeed(0.0f, 0.0f, 0.0f); + glm::vec3 apart(0.0f, 0.0f, 0.0f); + int rule1neighbors = 0; + int rule3neighbors = 0; + + for (int z = -1; z <= 1; z++) { + for (int y = -1; y <= 1; y++) { + + int neighborX = iX; + int neighborY = iY + y; + int neighborZ = iZ + z; + if (neighborX < 0 || neighborY < 0 || neighborZ < 0 || neighborX >= gridResolution || neighborY >= gridResolution || neighborZ >= gridResolution) { + continue; + } + + int gridIndex = gridIndex3Dto1D(neighborX, neighborY, neighborZ, gridResolution); + int startBoid = gridCellStartIndices[gridIndex]; + int endBoid = gridCellEndIndices[gridIndex]; + if (startBoid == -1 || endBoid == -1) { + continue; + } + + + for (int i = startBoid; i < endBoid; i++) { + + glm::vec3 neighborPos = pos[i]; + glm::vec3 neighborVel1 = vel1[i]; + + if (i == index) { + continue; + } + + float distance = glm::distance(neighborPos, pos[index]); + + if (distance < rule1Distance) { + averagePosition += neighborPos; + rule1neighbors++; + } + + + + if (distance < rule2Distance) { + apart -= (neighborPos - pos[index]); + } + + if (distance < rule3Distance) { + averageSpeed += neighborVel1; + rule3neighbors++; + } + + + + } + + + + } + } + + + if (rule1neighbors > 0) { + averagePosition /= static_cast(rule1neighbors); + velChange += (averagePosition - pos[index]) * rule1Scale; + } + + + velChange += (apart * rule2Scale); + + if (rule3neighbors > 0) { + averageSpeed /= static_cast(rule3neighbors); + velChange += (averageSpeed * rule3Scale); + } + + + glm::vec3 newVel = vel1[index] + velChange; + newVel = glm::clamp(newVel, maxSpeed * -1.0f, maxSpeed); + vel2[index] = newVel; + + + + + + + + } /** * Step the entire N-body simulation by `dt` seconds. */ void Boids::stepSimulationNaive(float dt) { - // TODO-1.2 - use the kernels you wrote to step the simulation forward in time. - // TODO-1.2 ping-pong the velocity buffers + // Use the kernels you wrote to step the simulation forward in time. + // ping-pong the velocity buffers + dim3 fullBlocksPerGrid((numObjects + blockSize - 1) / blockSize); + kernUpdateVelocityBruteForce<<>>(numObjects, dev_pos, dev_vel1, dev_vel2); + kernUpdatePos<<>>(numObjects, dt, dev_pos, dev_vel2); + glm::vec3* temp = dev_vel1; + dev_vel1 = dev_vel2; + dev_vel2 = temp; + } void Boids::stepSimulationScatteredGrid(float dt) { - // TODO-2.1 // Uniform Grid Neighbor search using Thrust sort. // In Parallel: // - label each particle with its array index as well as its grid index. @@ -374,11 +703,31 @@ void Boids::stepSimulationScatteredGrid(float dt) { // - Perform velocity updates using neighbor search // - Update positions // - Ping-pong buffers as needed + + dim3 fullBlocksPerGrid((numObjects + blockSize - 1) / blockSize); + + //Reset the start and end indices to -1 + kernResetIntBuffer << > > (gridCellCount, dev_gridCellStartIndices, -1); + kernResetIntBuffer << > > (gridCellCount, dev_gridCellEndIndices, -1); + + + //Sort based on grid indices + kernComputeIndices<<>>(numObjects, gridSideCount, gridMinimum, gridInverseCellWidth, dev_pos, dev_particleArrayIndices, dev_particleGridIndices); + thrust::sort_by_key(dev_thrust_particleGridIndices, dev_thrust_particleGridIndices + numObjects, dev_thrust_particleArrayIndices); + + kernIdentifyCellStartEnd << > > (numObjects, dev_particleGridIndices, dev_gridCellStartIndices, dev_gridCellEndIndices); + kernUpdateVelNeighborSearchScattered <<> > (numObjects, gridSideCount, gridMinimum, + gridInverseCellWidth, gridCellWidth, dev_gridCellStartIndices, dev_gridCellEndIndices, dev_particleArrayIndices, dev_pos, dev_vel1, dev_vel2); + kernUpdatePos << > > (numObjects, dt, dev_pos, dev_vel2); + + //ping pong buffers + glm::vec3* temp = dev_vel1; + dev_vel1 = dev_vel2; + dev_vel2 = temp; } void Boids::stepSimulationCoherentGrid(float dt) { - // TODO-2.3 - start by copying Boids::stepSimulationNaiveGrid - // Uniform Grid Neighbor search using Thrust sort on cell-coherent data. + // In Parallel: // - Label each particle with its array index as well as its grid index. // Use 2x width grids @@ -388,10 +737,41 @@ void Boids::stepSimulationCoherentGrid(float dt) { // cell's data pointers in the array of boid indices // - BIG DIFFERENCE: use the rearranged array index buffer to reshuffle all // the particle data in the simulation array. - // CONSIDER WHAT ADDITIONAL BUFFERS YOU NEED // - Perform velocity updates using neighbor search // - Update positions // - Ping-pong buffers as needed. THIS MAY BE DIFFERENT FROM BEFORE. + dim3 fullBlocksPerGrid((numObjects + blockSize - 1) / blockSize); + + //Reset the start and end indices to -1 + kernResetIntBuffer << > > (gridCellCount, dev_gridCellStartIndices, -1); + kernResetIntBuffer << > > (gridCellCount, dev_gridCellEndIndices, -1); + kernComputeIndices << > > (numObjects, gridSideCount, gridMinimum, gridInverseCellWidth, dev_pos, dev_particleArrayIndices, dev_particleGridIndices); + + //Sort based on grid indices + thrust::sort_by_key(dev_thrust_particleGridIndices, dev_thrust_particleGridIndices + numObjects, dev_thrust_particleArrayIndices); + kernIdentifyCellStartEnd << > > (numObjects, dev_particleGridIndices, dev_gridCellStartIndices, dev_gridCellEndIndices); + + //Reshuffle the pos and vel1 data into coherent order for better cache locality + kernReshuffleBoidData <<>> (numObjects, dev_particleArrayIndices, dev_pos, dev_vel1, dev_reshuffledPos, dev_reshuffledVel1); + + //Swap Pointers so that dev_pos and dev_vel1 are coherent + glm::vec3* tempPos = dev_pos; + dev_pos = dev_reshuffledPos; + dev_reshuffledPos = tempPos; + + glm::vec3* tempVel = dev_vel1; + dev_vel1 = dev_reshuffledVel1; + dev_reshuffledVel1 = tempVel; + + kernUpdateVelNeighborSearchCoherent <<> > (numObjects, gridSideCount, gridMinimum, + gridInverseCellWidth, gridCellWidth, dev_gridCellStartIndices, dev_gridCellEndIndices, dev_pos, dev_vel1, dev_vel2); + kernUpdatePos << > > (numObjects, dt, dev_pos, dev_vel2); + + //ping pong buffers + glm::vec3* temp = dev_vel1; + dev_vel1 = dev_vel2; + dev_vel2 = temp; + } void Boids::endSimulation() { @@ -399,7 +779,16 @@ void Boids::endSimulation() { cudaFree(dev_vel2); cudaFree(dev_pos); - // TODO-2.1 TODO-2.3 - Free any additional buffers here. + // additional buffers + cudaFree(dev_particleArrayIndices); + cudaFree(dev_particleGridIndices); + cudaFree(dev_gridCellStartIndices); + cudaFree(dev_gridCellEndIndices); + cudaFree(dev_reshuffledPos); + cudaFree(dev_reshuffledVel1); + + + } void Boids::unitTest() { diff --git a/src/main.cpp b/src/main.cpp index 9c917c0..8e35c06 100644 --- a/src/main.cpp +++ b/src/main.cpp @@ -23,11 +23,11 @@ // LOOK-2.1 LOOK-2.3 - toggles for UNIFORM_GRID and COHERENT_GRID #define VISUALIZE 1 -#define UNIFORM_GRID 0 -#define COHERENT_GRID 0 +#define UNIFORM_GRID 1 +#define COHERENT_GRID 1 // LOOK-1.2 - change this to adjust particle count in the simulation -const int N_FOR_VIS = 5000; +const int N_FOR_VIS = 1000000; const float DT = 0.2f; /**