Reaction Diffusion Simulation in GLSL
/* | |
* Copyright (c) 2016 David Wicks, sansumbrella.com | |
* All rights reserved. | |
* | |
* Redistribution and use in source and binary forms, with or | |
* without modification, are permitted provided that the following | |
* conditions are met: | |
* | |
* Redistributions of source code must retain the above copyright | |
* notice, this list of conditions and the following disclaimer. | |
* Redistributions in binary form must reproduce the above copyright | |
* notice, this list of conditions and the following disclaimer in the | |
* documentation and/or other materials provided with the distribution. | |
* | |
* THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS | |
* "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT | |
* LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR | |
* A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT | |
* HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, | |
* SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT | |
* LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, | |
* DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY | |
* THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT | |
* (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE | |
* OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. | |
*/ |
#version 330 core | |
/// | |
/// Some shared math utilities. | |
/// | |
float vmax(vec2 v) { | |
return max(v.x, v.y); | |
} | |
float vmax(vec3 v) { | |
return max(max(v.x, v.y), v.z); | |
} | |
/// reduce vector values to [0.0, 1.0] while retaining linear relationship | |
vec3 unbloom(vec3 v) { | |
float m = vmax(v); | |
if (m > 1.0) { | |
return v / m; | |
} | |
return v; | |
} | |
float sdf_square(vec2 center, float half_size, vec2 coord) { | |
float d = vmax(abs(center - coord)); | |
return step(d, half_size); | |
} | |
float sdf_circle(vec2 center, float radius, vec2 coord) { | |
float l = length(center - coord); | |
return step(l, radius); | |
} |
#include "cinder/app/App.h" | |
#include "cinder/app/RendererGl.h" | |
#include "cinder/gl/gl.h" | |
#include "cinder/Log.h" | |
#include "Watchdog.h" | |
#include "CinderImGui.h" | |
using namespace ci; | |
using namespace ci::app; | |
using namespace std; | |
/// | |
/// Demonstrates using Framebuffer Objects as data buffers to run a simulation using GLSL. | |
/// The two data buffers (source and destination) store floating point data and alternate being read from and written to. | |
/// The draw calls simply draw fullscreen quads to cover the buffers (in opengl that's just [-1, 1], [-1, 1] without any projection). | |
/// | |
/// This uses [WatchDog]( https://github.com/simongeilfus/Watchdog ) to automatically reload files. | |
/// If you don't have/want to download watchdog, you can remove the includes/usage of watchdog and | |
/// use key presses to reload shaders while the app is running. | |
/// | |
class PingPongShadersApp : public App { | |
public: | |
void setup() override; | |
void keyDown(KeyEvent event) override; | |
void update() override; | |
void draw() override; | |
void mouseDown(MouseEvent event) override; | |
void mouseUp(MouseEvent event) override; | |
void mouseDrag(MouseEvent event) override; | |
void loadSimulationShader(); | |
void loadRenderShader(); | |
void seedSimulation(); | |
void addSquare(); | |
private: | |
gl::FboRef _source, _destination; | |
gl::BatchRef _simulation_batch; | |
gl::BatchRef _render_batch; | |
vec2 _mouse_pos; | |
bool _drawing = false; | |
}; | |
void PingPongShadersApp::loadSimulationShader() | |
{ | |
try { | |
auto glsl = gl::GlslProg::create(loadAsset("shaders/simulate.vs"), loadAsset("shaders/simulate.fs")); | |
glsl->uniform("texture_0", 0); | |
glsl->uniform("pixel_size", vec2(1.0f) / vec2(_destination->getSize())); | |
_simulation_batch = gl::Batch::create(geom::Rect(Rectf(-1.0f, -1.0f, 1.0f, 1.0f)).texCoords(vec2(0, 0), vec2(1, 0), vec2(1, 1), vec2(0, 1)), glsl); | |
CI_LOG_I("Reloaded simulation shader."); | |
} | |
catch (const std::exception &exc) { | |
CI_LOG_W("Unable to load simulation shaders: " << exc.what()); | |
} | |
} | |
void PingPongShadersApp::loadRenderShader() | |
{ | |
try { | |
auto glsl = gl::GlslProg::create(loadAsset("shaders/render.vs"), loadAsset("shaders/render.fs")); | |
_render_batch = gl::Batch::create(geom::Rect(Rectf(-1.0f, -1.0f, 1.0f, 1.0f)).texCoords(vec2(0, 0), vec2(1, 0), vec2(1, 1), vec2(0, 1)), glsl); | |
CI_LOG_I("Reloaded render shader."); | |
} | |
catch (const std::exception &exc) { | |
CI_LOG_W("Unable to load render shaders: " << exc.what()); | |
} | |
} | |
void PingPongShadersApp::setup() | |
{ | |
ui::initialize(); | |
gl::Fbo::Format fmt; | |
// We only need 2 color channels to simulate chemicals A and B | |
// If your GPU doesn't like GL_RG32F, you can try some other values, like: GL_RGB32F, GL_RGB16F | |
fmt.setColorTextureFormat(gl::Texture::Format().internalFormat(GL_RG32F).wrap(GL_CLAMP_TO_EDGE)); | |
_source = gl::Fbo::create(512, 512, fmt); | |
_destination = gl::Fbo::create(512, 512, fmt); | |
seedSimulation(); | |
addSquare(); | |
wd::watch("shaders/simulate.*", [this] (const fs::path &path) { | |
loadSimulationShader(); | |
}); | |
wd::watch("shaders/render.*", [this] (const fs::path &path) { | |
loadRenderShader(); | |
}); | |
} | |
void PingPongShadersApp::seedSimulation() | |
{ | |
gl::ScopedFramebuffer fbo(_source); | |
gl::clear(Color(1, 0, 0)); // fill with chemical "a" | |
} | |
void PingPongShadersApp::addSquare() | |
{ | |
gl::ScopedFramebuffer fbo(_source); | |
gl::ScopedViewport viewport(ivec2(0), _source->getSize()); | |
gl::ScopedBlendAdditive blend; | |
gl::ScopedMatrices matrices; | |
gl::setMatricesWindow(_source->getSize()); | |
gl::ScopedColor color(Color(0, 1, 0)); // draw in chemical "b" | |
vec2 center = vec2(_source->getSize()) * 0.5f; | |
gl::drawSolidRect(Rectf(center - vec2(5), center + vec2(5))); | |
} | |
void PingPongShadersApp::mouseDown(cinder::app::MouseEvent event) | |
{ | |
_drawing = true; | |
_mouse_pos = event.getPos(); | |
} | |
void PingPongShadersApp::mouseUp(cinder::app::MouseEvent event) | |
{ | |
_drawing = false; | |
} | |
void PingPongShadersApp::mouseDrag(cinder::app::MouseEvent event) | |
{ | |
_mouse_pos = event.getPos(); | |
} | |
void PingPongShadersApp::keyDown(cinder::app::KeyEvent event) | |
{ | |
switch (event.getCode()) { | |
case KeyEvent::KEY_r: | |
loadSimulationShader(); | |
loadRenderShader(); | |
break; | |
case KeyEvent::KEY_c: | |
seedSimulation(); | |
break; | |
} | |
} | |
void PingPongShadersApp::update() | |
{ | |
if (_drawing) { | |
// add chemical b by drawing with the mouse | |
gl::ScopedFramebuffer fbo(_source); | |
gl::ScopedViewport viewport(ivec2(0), _source->getSize()); | |
gl::ScopedBlendAdditive blend; | |
gl::ScopedMatrices matrices; | |
gl::setMatricesWindow(getWindowSize()); | |
gl::ScopedColor color(Color(0, 1, 0)); // add in chemical "b" | |
gl::drawSolidCircle(_mouse_pos, 30.0f); | |
} | |
if (_simulation_batch) { | |
// Run the simulation many times each frame so we can watch it unfold faster. | |
for (auto i = 0; i < 12; i += 1) { | |
gl::ScopedFramebuffer fbo(_destination); | |
gl::clear(Color(0, 0, 0)); | |
gl::ScopedViewport viewport(ivec2(0), _destination->getSize()); | |
gl::ScopedTextureBind tex0(_source->getColorTexture(), 0); | |
_simulation_batch->draw(); | |
std::swap(_source, _destination); | |
} | |
} | |
} | |
void PingPongShadersApp::draw() | |
{ | |
gl::clear(Color(0, 0, 0)); | |
if (_render_batch) { | |
gl::ScopedTextureBind tex0(_source->getColorTexture(), 0); | |
_render_batch->draw(); | |
} | |
} | |
void prepare_settings(App::Settings *settings) | |
{ | |
settings->setWindowSize(960, 960); | |
} | |
CINDER_APP(PingPongShadersApp, RendererGl, prepare_settings) |
#version 330 core | |
#include "math.glsl" | |
uniform sampler2D texture_0; | |
uniform float ciElapsedSeconds; | |
in Vertex { | |
vec2 tex_coord; | |
} v; | |
out vec4 f_color; | |
void main() { | |
// Calculate gradient across surface | |
vec2 pixel_size = vec2(1.0 / 512.0); | |
float dbdx = (texture(texture_0, v.tex_coord + vec2(pixel_size.x, 0)).y - texture(texture_0, v.tex_coord + vec2(-pixel_size.x, 0)).y) / 2.0; | |
float dbdy = (texture(texture_0, v.tex_coord + vec2(0, pixel_size.y)).y - texture(texture_0, v.tex_coord + vec2(0, -pixel_size.y)).y) / 2.0; | |
float gradient_softness = 0.5 / pixel_size.x; | |
vec3 gx = normalize(vec3(pixel_size.x * gradient_softness, 0, dbdx)); | |
vec3 gy = normalize(vec3(0, pixel_size.y * gradient_softness, dbdy)); | |
vec3 normal = normalize(cross(gy, gx)); | |
float simulation_value = smoothstep(0.1, 0.2, texture(texture_0, v.tex_coord).y); | |
// float theta = ciElapsedSeconds * 0.5; | |
// vec3 light_direction = normalize(vec3(cos(theta), sin(theta), -1.0)); | |
vec3 light_direction = normalize(vec3(1.0, -0.33, -1.0)); | |
float diffuse = max(dot(normal, light_direction), 0.0); | |
vec3 reflection = reflect(light_direction, normal); | |
const vec3 eye_direction = normalize(vec3(0.0, 0.0, 1.0)); | |
float specular = max(dot(eye_direction, reflection), 0.0); | |
specular = pow(specular, 6.0); | |
// heightmap-shaded view | |
vec3 a = vec3(0.0, 0.0, 0.0); | |
vec3 b = vec3(1.1, 1.05, 1.0); | |
vec3 specular_color = vec3(0.33, 1.1, 1.4); | |
f_color = vec4(mix(a, b, diffuse) * simulation_value + specular_color * specular, 1.0); | |
// value-only black & white view | |
// f_color = vec4(vec3(simulation_value), 1); | |
// normal-map view | |
// f_color = vec4((normal + vec3(1.0)) / 2.0, 1.0); | |
// raw rgba (ab01) view. | |
// f_color = vec4(texture(texture_0, v.tex_coord).xy, 0, 1); | |
} |
#version 330 core | |
#include "math.glsl" | |
uniform sampler2D texture_0; | |
uniform vec2 pixel_size; | |
in Vertex { | |
vec2 tex_coord; | |
} v; | |
out vec4 f_color; | |
void main() { | |
vec2 uv = v.tex_coord; | |
// Read the simulation values at our current position. | |
vec2 ab = texture(texture_0, uv).xy; | |
// The laplace is a measure of the difference between a cell's value and the average of its neighbor values. | |
// We need to add everything up to calculate this. | |
vec2 laplace = vec2(0); | |
laplace += texture(texture_0, uv + pixel_size * vec2(-1, -1)).xy * 0.05; | |
laplace += texture(texture_0, uv + pixel_size * vec2(0, -1)).xy * 0.2; | |
laplace += texture(texture_0, uv + pixel_size * vec2(1, -1)).xy * 0.05; | |
laplace += texture(texture_0, uv + pixel_size * vec2(-1, 0)).xy * 0.2; | |
laplace += ab * (-1.0); | |
laplace += texture(texture_0, uv + pixel_size * vec2(1, 0)).xy * 0.2; | |
laplace += texture(texture_0, uv + pixel_size * vec2(-1, 1)).xy * 0.05; | |
laplace += texture(texture_0, uv + pixel_size * vec2(0, 1)).xy * 0.2; | |
laplace += texture(texture_0, uv + pixel_size * vec2(1, 1)).xy * 0.05; | |
float diffusion_a = mix(0.9, 1.03, uv.x); | |
float diffusion_b = mix(0.4, 0.58, uv.y); | |
float feed = mix(0.05, 0.06, uv.x); | |
float kill = mix(0.060, 0.064, uv.y); | |
float dt = 1.0; | |
float a = ab.x; | |
float b = ab.y; | |
// Calculate our new simulation values using the Gray-Scott equation, ignoring timestep. | |
float ap = clamp(a + ((diffusion_a * laplace.x) - (a * b * b) + (feed * (1.0 - a))) * dt, 0, 1); | |
float bp = clamp(b + ((diffusion_b * laplace.y) + (a * b * b) - ((kill + feed) * b)) * dt, 0, 1); | |
f_color = vec4(ap, bp, 0, 1); | |
// Optionally seed the simulation with an SDF shape here: | |
vec4 food_only = vec4(1, 0, 0, 1); | |
vec4 food_and_consumer = vec4(1, 1, 0, 1); | |
float square = sdf_square(vec2(0.5, 0.55), 0.1, v.tex_coord); | |
float circle = sdf_circle(vec2(0.6, 0.3), 0.2, v.tex_coord); | |
// f_color = mix(food_only, food_and_consumer, circle); | |
// f_color = mix(food_only, food_and_consumer, max(circle, square)); | |
} |
#version 330 core | |
in vec4 ciPosition; | |
in vec2 ciTexCoord0; | |
out Vertex { | |
vec2 tex_coord; | |
} v; | |
void main() { | |
gl_Position = ciPosition; | |
v.tex_coord = ciTexCoord0; | |
} |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
This comment has been minimized.
sansumbrella commentedNov 3, 2016
To run:
.cpp
source into yourMyNewApp.cpp
file.assets/
directory.Build and run.
This has only been tested on OSX, but should work cross-platform.