tadded findOverlaps in sphere.py (untested), do not run sim when testing for contacts - 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 4d7fc20cc3683d3021cc70b24449b4bb18c4fb31
(DIR) parent 83ab042edfff0746b4347f5e56942cfda9f39797
(HTM) Author: Anders Damsgaard <anders.damsgaard@geo.au.dk>
Date: Tue, 30 Sep 2014 15:07:54 +0200
added findOverlaps in sphere.py (untested), do not run sim when testing for
contacts
Diffstat:
M python/sphere.py | 42 ++++++++-----------------------
M src/main.cpp | 4 ++--
2 files changed, 13 insertions(+), 33 deletions(-)
---
(DIR) diff --git a/python/sphere.py b/python/sphere.py
t@@ -3658,38 +3658,18 @@ class sim:
def findOverlaps(self):
'''
- Find all particle-particle overlaps by a n^2 contact search. The
- particle pair indexes and the distance of the overlaps is saved in the
- object itself as the ``.pairs`` and ``.overlaps`` members.
- Warning: this function is EXTREMELY slow.
+ Find all particle-particle overlaps by a n^2 contact search, which is
+ done in C++. The particle pair indexes and the distance of the overlaps
+ is saved in the object itself as the ``.pairs`` and ``.overlaps``
+ members.
'''
-
- # Allocate big arrays instead of dynamically allocating during the
- # iterative process. Blank positions will be erased afterwards
- max_avg_contacts_per_particle = 16
- self.pairs = numpy.empty((self.np[0]*max_avg_contacts_per_particle, 2))
- self.overlaps = numpy.empty_like(self.pairs)
-
- # Loop through particle pairs
- p = 0
- for i in numpy.arange(self.np):
- print('%.2f%%, idx %d' % (float(i/self.np[0]), i))
- for j in numpy.arange(self.np):
-
- if i < j:
- x_ij = self.x[i,:] - self.x[j,:]
- x_ij_length = numpy.sqrt(x_ij.dot(x_ij))
-
- delta_n = x_ij_length - (self.radius[i] + self.radius[j])
-
- if delta_n < 0.0:
- self.pairs[p,:] = [i,j]
- self.overlaps[p] = delta_n
- p += 1
-
- # trim output arrays
- self.pairs = self.pairs[:p+1,:]
- self.overlaps = self.overlaps[:p+1]
+ self.writebin(verbose=False)
+ subprocess.call('cd .. && ./sphere --contacts input/' + self.sid
+ + '.bin > output/' + self.sid + '.contacts.txt', shell=True)
+ contactdata = numpy.loadtxt('output/' + self.sid + '.contacts.txt')
+ self.pairs = numpy.array((contactdata[:,0], contactdata[:,1]),
+ dtype=numpy.int32)
+ self.overlaps = numpy.array(contactdata[:,2])
def forcechains(self, lc=200.0, uc=650.0, outformat='png', disp='2d'):
(DIR) diff --git a/src/main.cpp b/src/main.cpp
t@@ -179,7 +179,7 @@ int main(const int argc, const char *argv[])
dem.printContacts();
// Render image if requested
- if (render == 1)
+ else if (render == 1)
dem.render(method, max_val, lower_cutoff);
// Otherwise, start iterating through time
t@@ -195,7 +195,7 @@ int main(const int argc, const char *argv[])
dem.printContacts();
// Render image if requested
- if (render == 1)
+ else if (render == 1)
dem.render(method, max_val, lower_cutoff);
// Otherwise, start iterating through time