diff options
-rw-r--r-- | channel_2d_gl_interop.py | 98 | ||||
-rw-r--r-- | channel_2d_streamlines_gl_interop.py | 6 | ||||
-rw-r--r-- | channel_3d_gl_interop.py | 10 | ||||
-rw-r--r-- | channel_3d_sdf_grid_fin_volumetric_rendering_gl_interop.py | 112 | ||||
-rw-r--r-- | channel_3d_volumetric_rendering_gl_interop.py | 112 | ||||
-rw-r--r-- | ldc_3d_gl_interop.py | 10 | ||||
-rw-r--r-- | shell.nix | 18 | ||||
-rw-r--r-- | simulation.py | 21 | ||||
-rw-r--r-- | symbolic/optimizations.py | 1 | ||||
-rw-r--r-- | template/kernel.mako | 2 | ||||
-rw-r--r-- | template/opengl.mako | 28 | ||||
-rw-r--r-- | template/particles.mako | 10 | ||||
-rw-r--r-- | template/sdf.cl.mako | 4 | ||||
-rw-r--r-- | template/streamline.mako | 10 | ||||
-rw-r--r-- | trugfeuer_2d_gl_interop.py | 52 | ||||
-rw-r--r-- | utility/opengl.py | 79 | ||||
-rw-r--r-- | utility/particles.py | 5 | ||||
-rw-r--r-- | utility/streamline.py | 5 |
18 files changed, 276 insertions, 307 deletions
diff --git a/channel_2d_gl_interop.py b/channel_2d_gl_interop.py index cae2f5b..8e52d78 100644 --- a/channel_2d_gl_interop.py +++ b/channel_2d_gl_interop.py @@ -2,7 +2,7 @@ import numpy from string import Template from simulation import Lattice, Geometry -from utility.opengl import MomentsVertexBuffer +from utility.opengl import MomentsTexture from utility.particles import Particles from symbolic.generator import LBM @@ -85,15 +85,32 @@ lbm = LBM(D2Q9) window = glut_window(fullscreen = False) -vertex_shader = shaders.compileShader(Template(""" +vertex_shader = shaders.compileShader(""" #version 430 -layout (location=0) in vec4 CellMoments; - -out vec3 color; +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), @@ -102,44 +119,34 @@ vec3 blueRedPalette(float x) { ); } -vec2 fluidVertexAtIndex(uint i) { - const float y = floor(float(i) / $size_x); - return vec2( - i - $size_x*y, - y - ); -} - -void main() { - const vec2 idx = fluidVertexAtIndex(gl_VertexID); - - gl_Position = projection * vec4( - idx.x, - idx.y, - 0., - 1. - ); - - if (CellMoments[3] > 0.0) { - color = blueRedPalette(CellMoments[3] / 0.2); +void main(){ + const vec2 sample_pos = unit(frag_pos); + const vec4 data = texture(moments, sample_pos); + if (data.w < 0.0) { + result.a = 1.0; + result.rgb = vec3(0.4); } else { - color = vec3(0.4,0.4,0.4); + result.a = 1.0; + result.rgb = blueRedPalette(data[3] / 0.2); } -}""").substitute({ - 'size_x': lattice_x, - 'inflow': inflow -}), GL_VERTEX_SHADER) +} +""").substitute({ + "size_x": lattice_x, + "size_y": lattice_y +}), GL_FRAGMENT_SHADER) -fragment_shader = shaders.compileShader(""" +particle_fragment_shader = shaders.compileShader(""" #version 430 in vec3 color; +layout(location = 0) out vec4 frag_color; + void main(){ - gl_FragColor = vec4(color.xyz, 0.0); + frag_color = vec4(color.xyz, 0.0); }""", GL_FRAGMENT_SHADER) -particle_shader = shaders.compileShader(Template(""" +particle_vertex_shader = shaders.compileShader(Template(""" #version 430 layout (location=0) in vec4 Particles; @@ -160,9 +167,11 @@ void main() { }""").substitute({}), GL_VERTEX_SHADER) shader_program = shaders.compileProgram(vertex_shader, fragment_shader) -particle_program = shaders.compileProgram(particle_shader, fragment_shader) projection_id = shaders.glGetUniformLocation(shader_program, 'projection') +particle_program = shaders.compileProgram(particle_vertex_shader, particle_fragment_shader) +particle_projection_id = shaders.glGetUniformLocation(particle_program, 'projection') + lattice = Lattice( descriptor = D2Q9, geometry = Geometry(lattice_x, lattice_y), @@ -176,11 +185,10 @@ lattice.apply_material_map( get_channel_material_map(lattice.geometry)) lattice.sync_material() -moments_vbo = MomentsVertexBuffer(lattice) +moments_texture = MomentsTexture(lattice) particles = Particles( lattice, - moments_vbo, numpy.mgrid[ lattice.geometry.size_x//20:2*lattice.geometry.size_x//20:100j, 4*lattice.geometry.size_y//9:5*lattice.geometry.size_y//9:100000/100j @@ -190,7 +198,8 @@ def on_display(): for i in range(0,updates_per_frame): lattice.evolve() - moments_vbo.collect() + lattice.update_moments() + moments_texture.collect() for i in range(0,updates_per_frame): particles.update(aging = False) @@ -200,13 +209,18 @@ def on_display(): glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT) shaders.glUseProgram(shader_program) - glUniformMatrix4fv(projection_id, 1, False, numpy.ascontiguousarray(projection)) - moments_vbo.bind() - glPointSize(point_size) - glDrawArrays(GL_POINTS, 0, lattice.geometry.volume) + glUniformMatrix4fv(projection_id, 1, False, numpy.asfortranarray(projection)) + moments_texture.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() shaders.glUseProgram(particle_program) - glUniformMatrix4fv(projection_id, 1, False, numpy.asfortranarray(projection)) + glUniformMatrix4fv(particle_projection_id, 1, False, numpy.asfortranarray(projection)) particles.bind() glPointSize(point_size) glDrawArrays(GL_POINTS, 0, particles.count) diff --git a/channel_2d_streamlines_gl_interop.py b/channel_2d_streamlines_gl_interop.py index c1aa1e8..544d40a 100644 --- a/channel_2d_streamlines_gl_interop.py +++ b/channel_2d_streamlines_gl_interop.py @@ -2,7 +2,6 @@ import numpy from string import Template from simulation import Lattice, Geometry -from utility.opengl import MomentsVertexBuffer from utility.streamline import Streamlines from symbolic.generator import LBM @@ -137,16 +136,15 @@ lattice.apply_material_map( get_channel_material_map(lattice.geometry)) lattice.sync_material() -moments_vbo = MomentsVertexBuffer(lattice) streamline_texture = Streamlines( - lattice, moments_vbo, + lattice, list(map(lambda y: [2, y*lattice.geometry.size_y//48], range(1,48)))) def on_display(): for i in range(0,updates_per_frame): lattice.evolve() - moments_vbo.collect() + lattice.update_moments() streamline_texture.update() glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT) diff --git a/channel_3d_gl_interop.py b/channel_3d_gl_interop.py index 535d595..3b93f86 100644 --- a/channel_3d_gl_interop.py +++ b/channel_3d_gl_interop.py @@ -17,7 +17,6 @@ from geometry.box import Box from geometry.cylinder import Cylinder from utility.projection import Projection, Rotation -from utility.opengl import MomentsVertexBuffer from utility.mouse import MouseDragMonitor, MouseScrollMonitor lattice_x = 256 @@ -135,8 +134,10 @@ fragment_shader = shaders.compileShader(""" in vec3 color; +layout(location = 0) out vec4 frag_color; + void main(){ - gl_FragColor = vec4(color.xyz, 0.0); + frag_color = vec4(color.xyz, 0.0); }""", GL_FRAGMENT_SHADER) lighting_vertex_shader = shaders.compileShader(""" @@ -211,11 +212,8 @@ primitives = list(map(lambda material: material[0], filter(lambda material: no lattice.apply_material_map(material_map) lattice.sync_material() -moments_vbo = MomentsVertexBuffer(lattice) - particles = Particles( lattice, - moments_vbo, numpy.mgrid[ 2*lattice.geometry.size_x//100:4*lattice.geometry.size_x//100:particle_count/10000j, lattice.geometry.size_y//16:15*lattice.geometry.size_y//16:100j, @@ -231,7 +229,7 @@ def on_display(): for i in range(0,updates_per_frame): lattice.evolve() - moments_vbo.collect() + lattice.update_moments() for i in range(0,updates_per_frame): particles.update(aging = True) diff --git a/channel_3d_sdf_grid_fin_volumetric_rendering_gl_interop.py b/channel_3d_sdf_grid_fin_volumetric_rendering_gl_interop.py index 8806029..13896a8 100644 --- a/channel_3d_sdf_grid_fin_volumetric_rendering_gl_interop.py +++ b/channel_3d_sdf_grid_fin_volumetric_rendering_gl_interop.py @@ -22,13 +22,13 @@ from utility.projection import Projection, Rotation from utility.opengl import MomentsTexture from utility.mouse import MouseDragMonitor, MouseScrollMonitor -lattice_x = 180 -lattice_y = 100 +lattice_x = 170 +lattice_y = 90 lattice_z = 100 updates_per_frame = 5 -inflow = 0.075 +inflow = 0.1 relaxation_time = 0.51 lbm = LBM(D3Q19) @@ -52,52 +52,66 @@ boundary = Template(""" ) grid_fin = """ -v = rotate_z(translate(v, v3(center.x/2, center.y, center.z)), -0.8); -float width = 1; -float angle = 0.64; - -return add( - sadd( - sub( - rounded(box(v, v3(5, 28, 38)), 1), - rounded(box(v, v3(6, 26, 36)), 1) - ), - cylinder(translate(v, v3(0,0,-45)), 5, 12), - 1 - ), - sintersect( - box(v, v3(5, 28, 38)), - add( - add( - box(rotate_x(v, angle), v3(10, width, 100)), - box(rotate_x(v, -angle), v3(10, width, 100)) +float sdf(vec3 v) { + v = rotate_z(translate(v, v3(center.x/2, center.y, center.z)), -0.6); + const float width = 1; + const float angle = 0.64; + + return add( + sadd( + sub( + rounded(box(v, v3(5, 28, 38)), 1), + rounded(box(v, v3(6, 26, 36)), 1) ), + cylinder(translate(v, v3(0,0,-45)), 5, 12), + 1 + ), + sintersect( + box(v, v3(5, 28, 38)), add( add( - add( - box(rotate_x(translate(v, v3(0,0,25)), angle), v3(10, width, 100)), - box(rotate_x(translate(v, v3(0,0,25)), -angle), v3(10, width, 100)) - ), - add( - box(rotate_x(translate(v, v3(0,0,-25)), angle), v3(10, width, 100)), - box(rotate_x(translate(v, v3(0,0,-25)), -angle), v3(10, width, 100)) - ) + box(rotate_x(v, angle), v3(10, width, 100)), + box(rotate_x(v, -angle), v3(10, width, 100)) ), add( add( - box(rotate_x(translate(v, v3(0,0,50)), angle), v3(10, width, 100)), - box(rotate_x(translate(v, v3(0,0,50)), -angle), v3(10, width, 100)) + add( + box(rotate_x(translate(v, v3(0,0,25)), angle), v3(10, width, 100)), + box(rotate_x(translate(v, v3(0,0,25)), -angle), v3(10, width, 100)) + ), + add( + box(rotate_x(translate(v, v3(0,0,-25)), angle), v3(10, width, 100)), + box(rotate_x(translate(v, v3(0,0,-25)), -angle), v3(10, width, 100)) + ) ), add( - box(rotate_x(translate(v, v3(0,0,-50)), angle), v3(10, width, 100)), - box(rotate_x(translate(v, v3(0,0,-50)), -angle), v3(10, width, 100)) + add( + box(rotate_x(translate(v, v3(0,0,50)), angle), v3(10, width, 100)), + box(rotate_x(translate(v, v3(0,0,50)), -angle), v3(10, width, 100)) + ), + add( + box(rotate_x(translate(v, v3(0,0,-50)), angle), v3(10, width, 100)), + box(rotate_x(translate(v, v3(0,0,-50)), -angle), v3(10, width, 100)) + ) ) ) - ) - ), - 2 - ) -); + ), + 2 + ) + ); +} + +float sdf_bounding(vec3 v) { + v = rotate_z(translate(v, v3(center.x/2, center.y, center.z)), -0.6); + const float width = 1; + const float angle = 0.64; + + return sadd( + rounded(box(v, v3(5, 28, 38)), 1), + cylinder(translate(v, v3(0,0,-45)), 5, 12), + 1 + ); +} """ def glut_window(fullscreen = False): @@ -133,8 +147,10 @@ fragment_shader = shaders.compileShader(""" in vec3 color; +layout(location = 0) out vec4 frag_color; + void main(){ - gl_FragColor = vec4(vec3(0.5), 1.0); + frag_color = vec4(vec3(0.5), 1.0); }""", GL_FRAGMENT_SHADER) raycast_vertex_shader = shaders.compileShader(""" @@ -173,11 +189,11 @@ const vec3 center = vec3(${size_x/2.5}, ${size_y/2}, ${size_z/2}); #define OBSTACLE_STEPS 16 vec3 v3(float x, float y, float z) { - return vec3(x,y,z); + return vec3(x,y,z); } vec2 v2(float x, float y) { - return vec2(x,y); + return vec2(x,y); } vec3 fabs(vec3 x) { @@ -190,9 +206,7 @@ float fabs(float x) { <%include file="template/sdf.lib.glsl.mako"/> -float sdf(vec3 v) { - ${sdf_source} -} +${sdf_source} vec3 sdf_normal(vec3 v) { return normalize(vec3( @@ -250,14 +264,14 @@ vec3 trace(vec3 pos, vec3 ray) { for (int i=0; i < RAYMARCH_STEPS; ++i) { const vec3 sample_pos = pos + i*delta*ray; const vec4 data = texture(moments, unit(sample_pos)); - if (sdf(sample_pos) > delta) { - color += 0.5/RAYMARCH_STEPS * palette(length(data.yzw) / ${inflow}); + if (sdf_bounding(sample_pos) > delta) { + color += 1.0/RAYMARCH_STEPS * palette(length(data.yzw) / ${inflow}); } else { const vec4 obstacle_color = trace_obstacle(sample_pos, ray, delta); if (obstacle_color.w == 1.0) { return color + obstacle_color.xyz; } else { - color += 0.5/RAYMARCH_STEPS * palette(length(data.yzw) / ${inflow}); + color += 1.0/RAYMARCH_STEPS * palette(length(data.yzw) / ${inflow}); } } } @@ -301,7 +315,7 @@ lattice = Lattice( lattice.setup_channel_with_sdf_obstacle(grid_fin) -moments_texture = MomentsTexture(lattice) +moments_texture = MomentsTexture(lattice, include_materials = False) projection = Projection(distance = 2*lattice_x) rotation = Rotation([-0.5*lattice_x, -0.5*lattice_y, -0.5*lattice_z]) diff --git a/channel_3d_volumetric_rendering_gl_interop.py b/channel_3d_volumetric_rendering_gl_interop.py index 1b149da..38d8f54 100644 --- a/channel_3d_volumetric_rendering_gl_interop.py +++ b/channel_3d_volumetric_rendering_gl_interop.py @@ -6,7 +6,6 @@ from mako.lookup import TemplateLookup from pathlib import Path from simulation import Lattice, Geometry -from utility.particles import Particles from symbolic.generator import LBM import symbolic.D3Q27 as D3Q27 @@ -26,7 +25,7 @@ lattice_x = 200 lattice_y = 64 lattice_z = 64 -updates_per_frame = 8 +updates_per_frame = 50 inflow = 0.01 relaxation_time = 0.51 @@ -52,34 +51,46 @@ boundary = Template(""" ) channel = """ -return add( - add( - ssub( - box(translate(v, v3(25,center.y,center.z)), v3(5,32,32)), - add( - sphere(translate(v, v3(20,0.5*center.y,1.5*center.z)), 10), - sphere(translate(v, v3(30,1.5*center.y,0.5*center.z)), 10) +float sdf(vec3 v) { + return add( + add( + ssub( + box(translate(v, v3(25,center.y,center.z)), v3(5,32,32)), + add( + sphere(translate(v, v3(20,0.5*center.y,1.5*center.z)), 10), + sphere(translate(v, v3(30,1.5*center.y,0.5*center.z)), 10) + ), + 2 ), - 2 + ssub( + box(translate(v, v3(85,center.y,center.z)), v3(5,32,32)), + add( + sphere(translate(v, v3(90,1.5*center.y,1.5*center.z)), 10), + sphere(translate(v, v3(80,0.5*center.y,0.5*center.z)), 10) + ), + 2 + ) ), ssub( - box(translate(v, v3(85,center.y,center.z)), v3(5,32,32)), + box(translate(v, v3(145,center.y,center.z)), v3(5,32,32)), add( - sphere(translate(v, v3(90,1.5*center.y,1.5*center.z)), 10), - sphere(translate(v, v3(80,0.5*center.y,0.5*center.z)), 10) + cylinder(rotate_y(translate(v, v3(145,1.5*center.y,0.5*center.z)), 1), 10, 10), + cylinder(rotate_y(translate(v, v3(145,0.5*center.y,1.5*center.z)), -1), 10, 10) ), 2 ) - ), - ssub( - box(translate(v, v3(145,center.y,center.z)), v3(5,32,32)), + ); +} + +float sdf_bounding(vec3 v) { + return add( add( - cylinder(rotate_y(translate(v, v3(145,1.5*center.y,0.5*center.z)), 1), 10, 10), - cylinder(rotate_y(translate(v, v3(145,0.5*center.y,1.5*center.z)), -1), 10, 10) + box(translate(v, v3(25,center.y,center.z)), v3(5,32,32)), + box(translate(v, v3(85,center.y,center.z)), v3(5,32,32)) ), - 2 - ) -); + box(translate(v, v3(145,center.y,center.z)), v3(5,32,32)) + ); +} """ def glut_window(fullscreen = False): @@ -115,8 +126,10 @@ fragment_shader = shaders.compileShader(""" in vec3 color; +layout(location = 0) out vec4 frag_color; + void main(){ - gl_FragColor = vec4(vec3(0.5), 1.0); + frag_color = vec4(vec3(0.5), 1.0); }""", GL_FRAGMENT_SHADER) raycast_vertex_shader = shaders.compileShader(""" @@ -136,6 +149,10 @@ void main() { raycast_fragment_shader = shaders.compileShader(Template(""" #version 430 +#define EPSILON 1e-1 +#define RAYMARCH_STEPS 64 +#define OBSTACLE_STEPS 16 + in vec3 frag_pos; uniform vec4 camera_pos; @@ -144,22 +161,15 @@ uniform sampler3D moments; out vec4 result; -vec3 unit(vec3 v) { - return vec3(v[0] / ${size_x}, v[1] / ${size_y}, v[2] / ${size_z}); -} - +const vec3 cuboid = vec3(${size_x}, ${size_y}, ${size_z}); const vec3 center = vec3(${size_x/2.5}, ${size_y/2}, ${size_z/2}); -#define EPSILON 1e-1 -#define RAYMARCH_STEPS 64 -#define OBSTACLE_STEPS 16 - vec3 v3(float x, float y, float z) { - return vec3(x,y,z); + return vec3(x,y,z); } vec2 v2(float x, float y) { - return vec2(x,y); + return vec2(x,y); } vec3 fabs(vec3 x) { @@ -172,9 +182,7 @@ float fabs(float x) { <%include file="template/sdf.lib.glsl.mako"/> -float sdf(vec3 v) { - ${sdf_source} -} +${sdf_source} vec3 sdf_normal(vec3 v) { return normalize(vec3( @@ -184,14 +192,23 @@ vec3 sdf_normal(vec3 v) { )); } +vec3 unit(vec3 v) { + return v / cuboid; +} + vec3 palette(float x) { return mix( - vec3(0.0, 0.0, 0.0), - vec3(1.0, 0.0, 0.0), + vec3(0.251, 0.498, 0.498), + vec3(0.502, 0.082, 0.082), x ); } +vec3 getVelocityColorAt(vec3 pos) { + const vec4 data = texture(moments, unit(pos)); + return palette(length(data.yzw) / ${2*inflow}); +} + float distanceToLattice(vec3 v) { return box(v, vec3(${size_x},${size_y},${size_z})); } @@ -201,7 +218,6 @@ float maxRayLength(vec3 origin, vec3 ray) { } vec4 trace_obstacle(vec3 origin, vec3 ray, float delta) { - vec3 color = vec3(0); vec3 sample_pos = origin; float ray_dist = 0.0; @@ -214,8 +230,9 @@ vec4 trace_obstacle(vec3 origin, vec3 ray, float delta) { } if (abs(sdf_dist) < EPSILON) { + const vec3 color = vec3(0.5); const vec3 n = normalize(sdf_normal(sample_pos)); - return vec4(color + abs(dot(n, ray)), 1.0); + return vec4(abs(dot(n, ray)) * color, 1.0); } else { sample_pos = origin + ray_dist*ray; } @@ -225,21 +242,24 @@ vec4 trace_obstacle(vec3 origin, vec3 ray, float delta) { } vec3 trace(vec3 pos, vec3 ray) { - const float delta = maxRayLength(pos, ray) / RAYMARCH_STEPS; + const float depth = maxRayLength(pos, ray); + const float delta = depth / RAYMARCH_STEPS; + const float gamma = 0.5 / RAYMARCH_STEPS; vec3 color = vec3(0.0); + vec3 sample_pos = pos; for (int i=0; i < RAYMARCH_STEPS; ++i) { - const vec3 sample_pos = pos + i*delta*ray; - const vec4 data = texture(moments, unit(sample_pos)); - if (sdf(sample_pos) > delta) { - color += 0.5/RAYMARCH_STEPS * palette(length(data.yzw) / ${inflow}); + sample_pos += delta*ray; + + if (sdf_bounding(sample_pos) > delta) { + color += gamma * getVelocityColorAt(sample_pos); } else { const vec4 obstacle_color = trace_obstacle(sample_pos, ray, delta); if (obstacle_color.w == 1.0) { return color + obstacle_color.xyz; } else { - color += 0.5/RAYMARCH_STEPS * palette(length(data.yzw) / ${inflow}); + color += gamma * getVelocityColorAt(sample_pos); } } } @@ -283,7 +303,7 @@ lattice = Lattice( lattice.setup_channel_with_sdf_obstacle(channel) -moments_texture = MomentsTexture(lattice) +moments_texture = MomentsTexture(lattice, include_materials = False) projection = Projection(distance = 2*lattice_x) rotation = Rotation([-0.5*lattice_x, -0.5*lattice_y, -0.5*lattice_z]) diff --git a/ldc_3d_gl_interop.py b/ldc_3d_gl_interop.py index a8f6f8d..575dc72 100644 --- a/ldc_3d_gl_interop.py +++ b/ldc_3d_gl_interop.py @@ -13,7 +13,6 @@ from OpenGL.GLUT import * from OpenGL.GL import shaders from utility.projection import Projection, Rotation -from utility.opengl import MomentsVertexBuffer from utility.mouse import MouseDragMonitor, MouseScrollMonitor lattice_x = 64 @@ -117,8 +116,10 @@ fragment_shader = shaders.compileShader(""" in vec3 color; +layout(location = 0) out vec4 frag_color; + void main(){ - gl_FragColor = vec4(color.xyz, 0.0); + frag_color = vec4(color.xyz, 0.0); }""", GL_FRAGMENT_SHADER) particle_program = shaders.compileProgram(particle_shader, fragment_shader) @@ -142,11 +143,8 @@ lattice.apply_material_map( get_cavity_material_map(lattice.geometry)) lattice.sync_material() -moments_vbo = MomentsVertexBuffer(lattice) - particles = Particles( lattice, - moments_vbo, numpy.mgrid[ 8*lattice.geometry.size_x//10:9*lattice.geometry.size_x//10:10j, lattice.geometry.size_y//10:9*lattice.geometry.size_y//10:particle_count/100j, @@ -162,7 +160,7 @@ def on_display(): for i in range(0,updates_per_frame): lattice.evolve() - moments_vbo.collect() + lattice.update_moments() for i in range(0,updates_per_frame): particles.update(aging = True) @@ -9,15 +9,27 @@ pkgs.stdenvNoCC.mkDerivation rec { packageOverrides = self: super: { pyopencl = super.pyopencl.overridePythonAttrs(old: rec { buildInputs = with pkgs; [ - opencl-headers ocl-icd python37Packages.pybind11 - libGLU_combined + opencl-headers ocl-icd + python3Packages.pybind11 + libGLU libGL ]; # Enable OpenGL integration and fix build preBuild = '' python configure.py --cl-enable-gl + echo "CL_PRETEND_VERSION = \"1.2\"" >> siteconf.py export HOME=/tmp/pyopencl ''; }); + + pyopengl = super.pyopengl.overridePythonAttrs(old: rec { + version = "3.1.0"; + src = self.fetchPypi { + pname = "PyOpenGL"; + inherit version; + sha256 = "1byxjj6a8rwzhxhjqlc588zdad2qwxdd7vlam2653ylll31waiwv"; + }; + }); + }; in pkgs.python3.override { inherit packageOverrides; }; @@ -32,7 +44,7 @@ pkgs.stdenvNoCC.mkDerivation rec { sha256 = "1p2459dqvgakywvy5d31818hix4kic6ks9j4m582ypxyk5wj1ksz"; }; - buildInputs = with pkgs.python37Packages; [ + buildInputs = with pkgs.python3Packages; [ numpy ]; diff --git a/simulation.py b/simulation.py index 535310f..e2174ee 100644 --- a/simulation.py +++ b/simulation.py @@ -153,8 +153,17 @@ class Lattice: self.platform = cl.get_platforms()[platform] if opengl: - self.context = cl.Context( - properties=[(cl.context_properties.PLATFORM, self.platform)] + get_gl_sharing_context_properties()) + try: + self.context = cl.Context( + properties = [ + (cl.context_properties.PLATFORM, self.platform) + ] + get_gl_sharing_context_properties()) + except: + self.context = cl.Context( + properties = [ + (cl.context_properties.PLATFORM, self.platform) + ] + get_gl_sharing_context_properties(), + devices = [ self.platform.get_devices()[0] ]) else: self.context = cl.Context( properties=[(cl.context_properties.PLATFORM, self.platform)]) @@ -252,9 +261,7 @@ class Lattice: def sync(self): self.queue.finish() - def get_moments(self): - moments = numpy.ndarray(shape=(self.descriptor.d+1, self.memory.volume), dtype=self.float_type[0]) - + def update_moments(self): if self.tick: self.program.collect_moments( self.queue, self.grid.size(), self.layout, self.memory.cl_pop_b, self.memory.cl_moments) @@ -262,6 +269,8 @@ class Lattice: self.program.collect_moments( self.queue, self.grid.size(), self.layout, self.memory.cl_pop_a, self.memory.cl_moments) + def get_moments(self): + moments = numpy.ndarray(shape=(self.descriptor.d+1, self.memory.volume), dtype=self.float_type[0]) + self.update_moments() cl.enqueue_copy(self.queue, moments, self.memory.cl_moments).wait(); - return moments diff --git a/symbolic/optimizations.py b/symbolic/optimizations.py index 93dad09..8e43cc8 100644 --- a/symbolic/optimizations.py +++ b/symbolic/optimizations.py @@ -1,6 +1,7 @@ from sympy import * from sympy.codegen.rewriting import ReplaceOptim +from sympy.simplify import cse_main expand_square = ReplaceOptim( lambda e: e.is_Pow and e.exp.is_integer and e.exp == 2, diff --git a/template/kernel.mako b/template/kernel.mako index dd8eaee..2e26124 100644 --- a/template/kernel.mako +++ b/template/kernel.mako @@ -105,4 +105,4 @@ __kernel void collect_moments(__global ${float_type}* f, % for i, expr in enumerate(moments_assignment): moments[${pop_offset(i)} + gid] = ${ccode(expr.rhs)}; % endfor -} +}
\ No newline at end of file diff --git a/template/opengl.mako b/template/opengl.mako index 47264e8..daf5f66 100644 --- a/template/opengl.mako +++ b/template/opengl.mako @@ -13,11 +13,22 @@ def moments_cell(): 2: '(int2)(get_global_id(0), get_global_id(1))', 3: '(int4)(get_global_id(0), get_global_id(1), get_global_id(2), 0)' }.get(descriptor.d) + +def neighbor_offset(c_i): + return { + 2: lambda: c_i[1]*memory.size_x + c_i[0], + 3: lambda: c_i[2]*memory.size_x*memory.size_y + c_i[1]*memory.size_x + c_i[0] + }.get(descriptor.d)() + %> -__kernel void collect_gl_moments(__global ${float_type}* f, - __global int* material, - __global float4* moments) +__kernel void collect_gl_moments_and_materials_to_texture(__global ${float_type}* f, + __global int* material, +% if descriptor.d == 2: + __write_only image2d_t moments) +% elif descriptor.d == 3: + __write_only image3d_t moments) +% endif { const unsigned int gid = ${gid()}; @@ -52,18 +63,9 @@ __kernel void collect_gl_moments(__global ${float_type}* f, data.w = -material[gid]; } - moments[gid] = data; + write_imagef(moments, ${moments_cell()}, data); } -<% -def neighbor_offset(c_i): - return { - 2: lambda: c_i[1]*memory.size_x + c_i[0], - 3: lambda: c_i[2]*memory.size_x*memory.size_y + c_i[1]*memory.size_x + c_i[0] - }.get(descriptor.d)() - -%> - __kernel void collect_gl_moments_to_texture(__global ${float_type}* f, % if descriptor.d == 2: __write_only image2d_t moments) diff --git a/template/particles.mako b/template/particles.mako index 48191d9..0b9a67e 100644 --- a/template/particles.mako +++ b/template/particles.mako @@ -1,4 +1,4 @@ -__kernel void update_particles(__global float4* moments, +__kernel void update_particles(__global float* moments, __global int* material, __global float4* particles, __global float4* init_particles, @@ -14,13 +14,11 @@ __kernel void update_particles(__global float4* moments, const unsigned int gid = floor(particle.z)*${memory.size_x*memory.size_y} + floor(particle.y)*${memory.size_x} + floor(particle.x); % endif - const float4 moment = moments[gid]; - if (material[gid] == 1 && particle.w < 1.0) { - particle.x += moment.y; - particle.y += moment.z; + particle.x += moments[${1*memory.volume}+gid]; + particle.y += moments[${2*memory.volume}+gid]; % if descriptor.d == 3: - particle.z += moment.w; + particle.z += moments[${3*memory.volume}+gid]; % endif particle.w += min(particle.x, particle.y) * aging; } else { diff --git a/template/sdf.cl.mako b/template/sdf.cl.mako index b98b35f..0a62f4b 100644 --- a/template/sdf.cl.mako +++ b/template/sdf.cl.mako @@ -13,9 +13,7 @@ __constant float3 center = (float3)(${geometry.size_x/2.5}, ${geometry.size_y/2} <%include file="sdf.lib.glsl.mako"/> -float sdf(vec3 v) { - ${sdf_src} -} +${sdf_src} __kernel void setup_channel_with_sdf_obstacle(__global int* material) { const unsigned x = get_global_id(0); diff --git a/template/streamline.mako b/template/streamline.mako index 73c80e2..adf5efd 100644 --- a/template/streamline.mako +++ b/template/streamline.mako @@ -31,7 +31,7 @@ float3 blueRedPalette(float x) { ); } -__kernel void draw_streamline(__global float4* moments, +__kernel void draw_streamline(__global float* moments, __global int* material, __global float2* origins, __read_write image2d_t streamlines) @@ -40,14 +40,12 @@ __kernel void draw_streamline(__global float4* moments, for (int i = 0; i < ${2*memory.size_x}; ++i) { const unsigned int gid = round(particle.y)*${memory.size_x} + round(particle.x); - const float4 moment = moments[gid]; - if (material[gid] != 1) { break; } - particle.x += 0.5 * moment.y / 0.01; - particle.y += 0.5 * moment.z / 0.01; + |