tresolve thread lock with contactmodel=2 - sphere - GPU-based 3D discrete element method algorithm with optional fluid coupling
(HTM) git clone git://src.adamsgaard.dk/sphere
(DIR) Log
(DIR) Files
(DIR) Refs
(DIR) LICENSE
---
(DIR) commit e0731b7097b1e026f653a9f6debdfd88f0b7f281
(DIR) parent 44093ce85fa78cecda72209b859fb2a3f3f610ed
(HTM) Author: Anders Damsgaard <anders@adamsgaard.dk>
Date: Wed, 15 Feb 2023 11:44:25 +0100
resolve thread lock with contactmodel=2
tthe cuda function __syncthreads() leads to undefined behavior when
invoked for diverged threads. The synchronization calls were originally
added to increase read performance, but this does not appear to be a
problem in current tests. This commit also removes debug output.
Diffstat:
M src/contactsearch.cuh | 7 -------
1 file changed, 0 insertions(+), 7 deletions(-)
---
(DIR) diff --git a/src/contactsearch.cuh b/src/contactsearch.cuh
t@@ -209,8 +209,6 @@ __device__ void findContactsInCell(
Float3 distmod = MAKE_FLOAT3(0.0f, 0.0f, 0.0f);
if (findDistMod(&targetCell, &distmod) != -1) {
- __syncthreads();
-
//// Check and process particle-particle collisions
// Calculate linear cell ID
t@@ -223,8 +221,6 @@ __device__ void findContactsInCell(
// Make sure cell is not empty
if (startIdx != 0xffffffff) {
- __syncthreads();
-
// Highest particle index in cell + 1
unsigned int endIdx = dev_cellEnd[cellID];
t@@ -353,15 +349,12 @@ __global__ void topology(
gridPos.y = floor((x_a.y - devC_grid.origo[1]) / (devC_grid.L[1]/devC_grid.num[1]));
gridPos.z = floor((x_a.z - devC_grid.origo[2]) / (devC_grid.L[2]/devC_grid.num[2]));
- printf("\ntopology: idx_a = %d, gridpos = [%d, %d, %d]\n", idx_a, gridPos.x, gridPos.y, gridPos.z);
-
// Find overlaps between particle no. idx and all particles
// from its own cell + 26 neighbor cells
for (int z_dim=-1; z_dim<2; ++z_dim) { // z-axis
for (int y_dim=-1; y_dim<2; ++y_dim) { // y-axis
for (int x_dim=-1; x_dim<2; ++x_dim) { // x-axis
targetPos = gridPos + make_int3(x_dim, y_dim, z_dim);
- printf("\n targetPos = [%d, %d, %d]\n", targetPos.x, targetPos.y, targetPos.z);
findContactsInCell(targetPos, idx_a, x_a, radius_a,
dev_x_sorted,
dev_cellStart, dev_cellEnd,