tadded 'color' array for color index - 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 292006e5d71f0d6cf10d53580aacf121dd0affab
(DIR) parent 3126bc3283216bf1288e04d95a53973c384d09f4
(HTM) Author: Anders Damsgaard <anders.damsgaard@geo.au.dk>
Date: Thu, 19 Jun 2014 10:26:58 +0200
added 'color' array for color index
Diffstat:
M python/sphere.py | 27 +++++++++++++++++++++++----
M src/constants.h | 2 +-
M src/datatypes.h | 1 +
M src/file_io.cpp | 10 ++++++++--
M src/sphere.cpp | 1 +
5 files changed, 34 insertions(+), 7 deletions(-)
---
(DIR) diff --git a/python/sphere.py b/python/sphere.py
t@@ -12,7 +12,7 @@ numpy.seterr(all='warn', over='raise')
# Sphere version number. This field should correspond to the value in
# `../src/constants.h`.
-VERSION=1.01
+VERSION=1.02
class sim:
'''
t@@ -324,6 +324,8 @@ class sim:
# The number of DEM time steps to perform between CFD updates
self.ndem = numpy.array(1)
+ # Particle color marker
+ self.color = numpy.zeros(self.np, dtype=numpy.int32)
def __cmp__(self, other):
'''
t@@ -580,6 +582,10 @@ class sim:
print(83)
return 83
+ if ((self.color != other.color)).any():
+ print(90)
+ return 90
+
# All equal
return 0
t@@ -597,7 +603,8 @@ class sim:
es = numpy.zeros(1),
ev_dot = numpy.zeros(1),
ev = numpy.zeros(1),
- p = numpy.zeros(1)):
+ p = numpy.zeros(1),
+ color = 0):
'''
Add a single particle to the simulation object. The only required
parameters are the position (x) and the radius (radius).
t@@ -646,6 +653,7 @@ class sim:
self.ev_dot = numpy.append(self.ev_dot, ev_dot)
self.ev = numpy.append(self.ev, ev)
self.p = numpy.append(self.p, p)
+ self.color = numpy.append(self.color, color)
def deleteAllParticles(self):
'''
t@@ -666,6 +674,7 @@ class sim:
self.ev_dot = numpy.zeros(self.np, dtype=numpy.float64)
self.ev = numpy.zeros(self.np, dtype=numpy.float64)
self.p = numpy.zeros(self.np, dtype=numpy.float64)
+ self.color = numpy.zeros(self.np, dtype=numpy.int32)
def readbin(self, targetbin, verbose = True, bonds = True, devsmod = True,
esysparticle = False):
t@@ -907,6 +916,14 @@ class sim:
self.maxiter = numpy.fromfile(fh, dtype=numpy.uint32, count=1)
if (self.version >= 1.01):
self.ndem = numpy.fromfile(fh, dtype=numpy.uint32, count=1)
+ else:
+ self.ndem = 1
+
+ if (self.version >= 1.02):
+ self.color =\
+ numpy.fromfile(fh, dtype=numpy.int32, count=self.np)
+ else:
+ self.color = numpy.zeros(self.np, dtype=numpy.int32)
finally:
if fh is not None:
t@@ -1054,6 +1071,8 @@ class sim:
fh.write(self.maxiter.astype(numpy.uint32))
fh.write(self.ndem.astype(numpy.uint32))
+ fh.write(self.color.astype(numpy.int32))
+
finally:
if fh is not None:
fh.close()
t@@ -2259,6 +2278,7 @@ class sim:
self.angvel[I,2] = 0.0
self.vel[I,0] = 0.0 # x-dim
self.vel[I,1] = 0.0 # y-dim
+ self.color[I] = -1
# Fix horizontal velocity to specific value of uppermost particles
d_max_top = numpy.max(self.radius[numpy.nonzero(self.x[:,2] >
t@@ -2270,6 +2290,7 @@ class sim:
self.angvel[I,2] = 0.0
self.vel[I,0] = (z_max-z_min)*shear_strain_rate
self.vel[I,1] = 0.0 # y-dim
+ self.color[I] = -1
# Set wall tangential viscosity to zero
self.gamma_wt[0] = 0.0
t@@ -2304,7 +2325,6 @@ class sim:
:type total: float
'''
-
# Computational time step (O'Sullivan et al, 2003)
#self.time_dt[0] = 0.17 * \
#math.sqrt((4.0/3.0 * math.pi * r_min**3 * self.rho[0]) \
t@@ -2571,7 +2591,6 @@ class sim:
# Debonding distance
self.db[0] = (1.0 + theta/2.0) * self.V_b**(1.0/3.0)
-
def bond(self, i, j):
'''
Create a bond between particles with index i and j
(DIR) diff --git a/src/constants.h b/src/constants.h
t@@ -16,7 +16,7 @@ const Float PI = 3.14159265358979;
const unsigned int ND = 3;
// Define source code version
-const double VERSION = 1.01;
+const double VERSION = 1.02;
// Max. number of contacts per particle
//const int NC = 16;
(DIR) diff --git a/src/datatypes.h b/src/datatypes.h
t@@ -29,6 +29,7 @@ struct Kinematics {
uint2 *bonds; // Particle bond pairs
Float4 *bonds_delta; // Particle bond displacement
Float4 *bonds_omega; // Particle bond rotation
+ int *color; // Color index for visualization
};
// Structure containing individual particle energies
(DIR) diff --git a/src/file_io.cpp b/src/file_io.cpp
t@@ -8,6 +8,7 @@
#include "constants.h"
#include "sphere.h"
+
// Get the address of the first byte of an object's representation
// See Stroustrup (2008) p. 388
template<class T>
t@@ -91,6 +92,7 @@ void DEM::readbin(const char *target)
k.angpos = new Float4[np];
k.angvel = new Float4[np];
k.torque = new Float4[np];
+ k.color = new int[np];
e.es_dot = new Float[np];
e.es = new Float[np];
t@@ -295,11 +297,12 @@ void DEM::readbin(const char *target)
cout << "Done" << std::endl;
}
+ for (i = 0; i<np; ++i)
+ ifs.read(as_bytes(k.color[i]), sizeof(int));
+
// Close file if it is still open
if (ifs.is_open())
ifs.close();
-
-
}
// Write DEM data to binary file
t@@ -491,6 +494,9 @@ void DEM::writebin(const char *target)
ofs.write(as_bytes(ns.ndem), sizeof(unsigned int));
}
+ for (i = 0; i<np; ++i)
+ ofs.write(as_bytes(k.color[i]), sizeof(int));
+
// Close file if it is still open
if (ofs.is_open())
ofs.close();
(DIR) diff --git a/src/sphere.cpp b/src/sphere.cpp
t@@ -88,6 +88,7 @@ DEM::~DEM(void)
delete[] k.angpos;
delete[] k.angvel;
delete[] k.torque;
+ delete[] k.color;
delete[] e.es_dot;
delete[] e.es;
delete[] e.ev_dot;