diff options
Diffstat (limited to 'lid_driven_cavity/opencl_gl_interop/ldc_2d.py')
-rw-r--r-- | lid_driven_cavity/opencl_gl_interop/ldc_2d.py | 180 |
1 files changed, 180 insertions, 0 deletions
diff --git a/lid_driven_cavity/opencl_gl_interop/ldc_2d.py b/lid_driven_cavity/opencl_gl_interop/ldc_2d.py new file mode 100644 index 0000000..daac486 --- /dev/null +++ b/lid_driven_cavity/opencl_gl_interop/ldc_2d.py @@ -0,0 +1,180 @@ +import numpy +import time +from string import Template + +from boltzgen import Generator, Geometry +from boltzgen.lbm.lattice import D2Q9 +from boltzgen.lbm.model import BGK + +from common import CellList + +from OpenGL.GL import * +from OpenGL.GLUT import * +from OpenGL.GL import shaders +from pyrr import matrix44 + +geometry = Geometry(512, 512) + +functions = ['collide_and_stream', 'equilibrilize', 'collect_moments', 'momenta_boundary'] +extras = ['cell_list_dispatch', 'opencl_gl_interop'] + +precision = 'single' +streaming = 'AA' + +import AA +import AB + +def glut_window(fullscreen = False): + glutInit(sys.argv) + glutInitDisplayMode(GLUT_RGBA | GLUT_DOUBLE | GLUT_DEPTH) + + if fullscreen: + window = glutEnterGameMode() + else: + glutInitWindowSize(800, 500) + glutInitWindowPosition(0, 0) + window = glutCreateWindow("LDC 2D") + + return window + +window = glut_window(fullscreen = False) + +Lattice = eval('%s.Lattice' % streaming) +HelperTemplate = eval('%s.HelperTemplate' % streaming) +MomentsTexture = eval('%s.MomentsTexture' % streaming) + +generator = Generator( + model = BGK(D2Q9, tau = 0.53), + target = 'cl', + precision = precision, + streaming = streaming, + index = 'ZYX', + layout = 'SOA') + +kernel_src = generator.kernel(geometry, functions, extras) +kernel_src += generator.custom(geometry, HelperTemplate) + +lattice = Lattice(geometry, kernel_src, D2Q9, precision = precision) +moments = MomentsTexture(lattice) + +gid = lattice.memory.gid + +ghost_cells = CellList(lattice.context, lattice.queue, lattice.float_type, + [ gid(x,y) for x, y in geometry.cells() if x == 0 or y == 0 or x == geometry.size_x-1 or y == geometry.size_y-1 ]) +bulk_cells = CellList(lattice.context, lattice.queue, lattice.float_type, + [ gid(x,y) for x, y in geometry.inner_cells() if x > 1 and x < geometry.size_x-2 and y > 1 and y < geometry.size_y-2 ]) +wall_cells = CellList(lattice.context, lattice.queue, lattice.float_type, + [ gid(x,y) for x, y in geometry.inner_cells() if x == 1 or y == 1 or x == geometry.size_x-2 ]) +lid_cells = CellList(lattice.context, lattice.queue, lattice.float_type, + [ gid(x,y) for x, y in geometry.inner_cells() if y == geometry.size_y-2 ]) + +if streaming == 'AB': + lattice.schedule('collide_and_stream_cells', bulk_cells) + lattice.schedule('velocity_momenta_boundary_cells', wall_cells, numpy.array([0.0, 0.0], dtype=lattice.float_type[0])) + lattice.schedule('velocity_momenta_boundary_cells', lid_cells, numpy.array([0.1, 0.0], dtype=lattice.float_type[0])) + +elif streaming == 'AA': + lattice.schedule_tick('collide_and_stream_cells_tick', bulk_cells) + lattice.schedule_tick('velocity_momenta_boundary_cells_tick', wall_cells, numpy.array([0.0, 0.0], dtype=lattice.float_type[0])) + lattice.schedule_tick('velocity_momenta_boundary_cells_tick', lid_cells, numpy.array([0.1, 0.0], dtype=lattice.float_type[0])) + + lattice.schedule_tock('equilibrilize_cells_tick', ghost_cells) + lattice.schedule_tock('collide_and_stream_cells_tock', bulk_cells) + lattice.schedule_tock('velocity_momenta_boundary_cells_tock', wall_cells, numpy.array([0.0, 0.0], dtype=lattice.float_type[0])) + lattice.schedule_tock('velocity_momenta_boundary_cells_tock', lid_cells, numpy.array([0.1, 0.0], dtype=lattice.float_type[0])) + + +def get_projection(width, height): + world_height = geometry.size_y + world_width = world_height / height * width + + projection = matrix44.create_orthogonal_projection(-world_width/2, world_width/2, -world_height/2, world_height/2, -1, 1) + translation = matrix44.create_from_translation([-geometry.size_x/2, -geometry.size_y/2, 0]) + + point_size = width / world_width + + return numpy.matmul(translation, projection), point_size + +vertex_shader = shaders.compileShader(""" +#version 430 + +layout (location=0) in vec4 vertex; + out vec2 frag_pos; + +uniform mat4 projection; + +void main() { + gl_Position = projection * vertex; + frag_pos = vertex.xy; +}""", GL_VERTEX_SHADER) + +fragment_shader = shaders.compileShader(Template(""" +#version 430 + +in vec2 frag_pos; + +uniform sampler2D moments; + +out vec4 result; + +vec2 unit(vec2 v) { + return vec2(v[0] / $size_x, v[1] / $size_y); +} + +vec3 blueRedPalette(float x) { + return mix( + vec3(0.0, 0.0, 1.0), + vec3(1.0, 0.0, 0.0), + x + ); +} + +void main(){ + const vec2 sample_pos = unit(frag_pos); + const vec4 data = texture(moments, sample_pos); + result.a = 1.0; + result.rgb = blueRedPalette(data[3] / 0.1); +} +""").substitute({ + "size_x": geometry.size_x, + "size_y": geometry.size_y, +}), GL_FRAGMENT_SHADER) + +shader_program = shaders.compileProgram(vertex_shader, fragment_shader) +projection_id = shaders.glGetUniformLocation(shader_program, 'projection') + +def on_display(): + for i in range(0,100): + lattice.evolve() + + moments.collect() + + glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT) + + shaders.glUseProgram(shader_program) + glUniformMatrix4fv(projection_id, 1, False, numpy.asfortranarray(projection)) + moments.bind() + + glBegin(GL_POLYGON) + glVertex(0,0,0) + glVertex(lattice.geometry.size_x,0,0) + glVertex(lattice.geometry.size_x,lattice.geometry.size_y,0) + glVertex(0,lattice.geometry.size_y,0) + glEnd() + + glutSwapBuffers() + +def on_reshape(width, height): + global projection, point_size + glViewport(0,0,width,height) + projection, point_size = get_projection(width, height) + +def on_timer(t): + glutTimerFunc(t, on_timer, t) + glutPostRedisplay() + +glutDisplayFunc(on_display) +glutReshapeFunc(on_reshape) +glutTimerFunc(10, on_timer, 10) + +glutMainLoop() |