aboutsummaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--implosion.py28
1 files changed, 20 insertions, 8 deletions
diff --git a/implosion.py b/implosion.py
index 25d53c2..8a3e530 100644
--- a/implosion.py
+++ b/implosion.py
@@ -68,7 +68,8 @@ float equilibrium(float d, float2 v, int i, int j) {
}
__kernel void collide_and_stream(__global float* pop_a,
- __global float* pop_b)
+ __global float* pop_b,
+ __global int* material)
{
const unsigned int gid = get_global_id(0);
const uint2 cell = cellAtGid(gid);
@@ -76,7 +77,7 @@ __kernel void collide_and_stream(__global float* pop_a,
float d = density(&pop_b[gid*9]);
float2 v = velocity(&pop_b[gid*9],d);
- if ( cell.x >= 2 && cell.x <= $nX-3 && cell.y >= 2 && cell.y <= $nY-3 ) {
+ if ( material[gid] == 1 ) {
for ( int i = -1; i <= 1; ++i ) {
for ( int j = -1; j <= 1; ++j ) {
pop_a[gidOfCell(cell.x+i, cell.y+j)*9 + indexOfDirection(i,j)] =
@@ -84,9 +85,7 @@ __kernel void collide_and_stream(__global float* pop_a,
}
}
}
- else if ( ((cell.y == 1 || cell.y == $nY-2) && (cell.x > 0 && cell.x < $nX-1)) ||
- ((cell.x == 1 || cell.x == $nX-2) && (cell.y > 0 && cell.y < $nY-1)) )
- {
+ else if ( material[gid] == 2 ) {
for ( int i = -1; i <= 1; ++i ) {
for ( int j = -1; j <= 1; ++j ) {
pop_a[gidOfCell(cell.x-i, cell.y-j)*9 + indexOfDirection(-i,-j)] =
@@ -112,15 +111,28 @@ class D2Q9_BGK_Lattice:
self.np_pop_a = numpy.ndarray(shape=(self.nCells, 9), dtype=numpy.float32)
self.np_pop_b = numpy.ndarray(shape=(self.nCells, 9), dtype=numpy.float32)
+ self.np_material = numpy.ndarray(shape=(self.nCells, 1), dtype=numpy.int32)
+
+ self.setup_geometry()
self.equilibrilize()
self.setup_anomaly()
self.cl_pop_a = cl.Buffer(self.context, mf.READ_WRITE | mf.USE_HOST_PTR, hostbuf=self.np_pop_a)
self.cl_pop_b = cl.Buffer(self.context, mf.READ_WRITE | mf.USE_HOST_PTR, hostbuf=self.np_pop_b)
+ self.cl_material = cl.Buffer(self.context, mf.READ_ONLY | mf.USE_HOST_PTR, hostbuf=self.np_material)
self.build_kernel()
+ def setup_geometry(self):
+ self.np_material[:] = 0
+ for x in range(1,self.nX-1):
+ for y in range(1,self.nY-1):
+ if x == 1 or y == 1 or x == self.nX-2 or y == self.nY-2:
+ self.np_material[self.idx(x,y)] = 2
+ else:
+ self.np_material[self.idx(x,y)] = 1
+
def equilibrilize(self):
self.np_pop_a[:,:] = [ 1./36., 1./9., 1./36., 1./9. , 4./9., 1./9. , 1./36 , 1./9., 1./36. ]
self.np_pop_b[:,:] = [ 1./36., 1./9., 1./36., 1./9. , 4./9., 1./9. , 1./36 , 1./9., 1./36. ]
@@ -148,11 +160,11 @@ class D2Q9_BGK_Lattice:
def evolve(self):
if self.tick:
self.tick = False
- self.program.collide_and_stream(self.queue, (self.nCells,), None, self.cl_pop_a, self.cl_pop_b)
+ self.program.collide_and_stream(self.queue, (self.nCells,), None, self.cl_pop_a, self.cl_pop_b, self.cl_material)
self.queue.finish()
else:
self.tick = True
- self.program.collide_and_stream(self.queue, (self.nCells,), None, self.cl_pop_b, self.cl_pop_a)
+ self.program.collide_and_stream(self.queue, (self.nCells,), None, self.cl_pop_b, self.cl_pop_a, self.cl_material)
self.queue.finish()
def show(self, i):
@@ -180,7 +192,7 @@ nUpdates = 100
start = timer()
-for i in range(1,nUpdates):
+for i in range(0,nUpdates):
LBM.evolve()
end = timer()