lisyarus/psemek
1
0
Fork 1
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

More 2d shallow water tests

Browse source
This commit is contained in:
Nikita Lisitsa 2025-01-30 13:07:30 +03:00
parent fa9c5efaf0
commit a6268663da
1 changed files with 114 additions and 22 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

136
examples/water_2d.cpp
View file

@ -2,6 +2,7 @@
#include <psemek/app/default_application_factory.hpp> #include <psemek/app/default_application_factory.hpp>
#include <psemek/gfx/painter.hpp> #include <psemek/gfx/painter.hpp>
#include <psemek/math/orthographic.hpp> #include <psemek/math/orthographic.hpp>
#include <psemek/math/camera.hpp>
#include <psemek/log/log.hpp> #include <psemek/log/log.hpp>
#include <psemek/util/array.hpp> #include <psemek/util/array.hpp>
#include <psemek/random/device.hpp> #include <psemek/random/device.hpp>
@ -31,7 +32,7 @@ struct water_2d_app
random::generator rng{random::device{}}; random::generator rng{random::device{}};
int const N = 128; int const N = 256;
float time = 0.f; float time = 0.f;
@ -48,14 +49,23 @@ struct water_2d_app
bool paused = true; bool paused = true;
bool show_water = true; bool show_water = true;
bool show_velocity = false; bool show_velocity = false;
bool show_particles = false;
bool erosion_on = false; bool erosion_on = false;
bool rain_on = false; bool rain_on = false;
struct particle
{
math::point<float, 2> position;
int lifetime;
};
std::vector<particle> particles;
}; };
water_2d_app::water_2d_app(options const &, context const & ctx) water_2d_app::water_2d_app(options const &, context const & ctx)
{ {
(void)ctx; (void)ctx;
ctx.vsync(true); // ctx.vsync(true);
simulation_area[0] = {-1.f, 1.f}; simulation_area[0] = {-1.f, 1.f};
simulation_area[1] = {-1.f, 1.f}; simulation_area[1] = {-1.f, 1.f};
@ -86,8 +96,8 @@ water_2d_app::water_2d_app(options const &, context const & ctx)
// if (x == N / 2 && (y != N / 4 && y != 3 * N / 4)) // if (x == N / 2 && (y != N / 4 && y != 3 * N / 4))
// bed(x, y) = 100.f; // bed(x, y) = 100.f;
float tx = (x + 0.5f) / N; [[maybe_unused]] float tx = (x + 0.5f) / N;
float ty = (y + 0.5f) / N; [[maybe_unused]] float ty = (y + 0.5f) / N;
float n = 0.f; float n = 0.f;
@ -106,7 +116,9 @@ water_2d_app::water_2d_app(options const &, context const & ctx)
// n = 4.f * std::abs(2.f * noise(tx, ty) - 1.f) * std::pow(4.f * tx * (1.f - tx), 2.f); // n = 4.f * std::abs(2.f * noise(tx, ty) - 1.f) * std::pow(4.f * tx * (1.f - tx), 2.f);
// River // River
n = std::min(4.f, math::lerp(10.f, 20.f, noise(tx, ty)) * std::abs(ty - math::lerp({0.4f, 0.6f}, noise(tx, 0.f)))); // n = std::min(4.f, math::lerp(10.f, 20.f, noise(tx, ty)) * std::abs(ty - math::lerp({0.4f, 0.6f}, noise(tx, 0.f))));
// n = 0.f;
// Canyon // Canyon
// n = 10.f * math::clamp(10.f * math::lerp(-1.f, 1.f, noise(tx, ty)), {0.f, 1.f}); // n = 10.f * math::clamp(10.f * math::lerp(-1.f, 1.f, noise(tx, ty)), {0.f, 1.f});
@ -128,13 +140,17 @@ water_2d_app::water_2d_app(options const &, context const & ctx)
// if (x == 0) // if (x == 0)
// water(x, y) = 100.f - bed(x, y); // water(x, y) = 100.f - bed(x, y);
water(x, y) = std::max(0.f, 1.f - bed(x, y)) * (0.5f + 0.5f * noise(tx, ty));
flowx(x, y) = (2.f * ty - 1.f) * water(x, y) * 10.f / N;
} }
} }
// water(N / 2, N / 2) = 10000.f; // water(N / 2, N / 2) = 10000.f;
} }
float const dt = 0.002f; float const dt = 0.001f;
void water_2d_app::update() void water_2d_app::update()
{ {
@ -144,11 +160,15 @@ void water_2d_app::update()
float const dx = simulation_area[0].length() / N; float const dx = simulation_area[0].length() / N;
float const dy = simulation_area[1].length() / N; float const dy = simulation_area[1].length() / N;
float const g = 10.f; float const g = 10.f;
float const friction = std::pow(0.5f, dt); float const friction = std::pow(1.f, dt);
float const viscosity = 0.f; float const viscosity = 0.f;
float const sediment_capacity = 0.1f; float const sediment_capacity = 0.1f;
float const erosion = 1.f; float const erosion = 1.f;
float const deposition = 10.f; float const deposition = 10.f;
float const coriolis = 0.01f;
int const max_particles = 16 * 1024;
int const particles_per_frame = 16;
int const max_particle_lifetime = max_particles / particles_per_frame;
// Init boundary flows // Init boundary flows
for (int x = 0; x < N; ++x) for (int x = 0; x < N; ++x)
@ -160,13 +180,16 @@ void water_2d_app::update()
// if (x < N / 4) flowy(x, N) = -3.f / N; // if (x < N / 4) flowy(x, N) = -3.f / N;
// flowy(x, 0) = 5.f * std::sin(10.f * time) / N; // flowy(x, 0) = 5.f * std::sin(10.f * time) / N;
flowy(x, 0) = 0.f;
flowy(x, N) = 0.f;
} }
for (int y = 0; y < N; ++y) for (int y = 0; y < N; ++y)
{ {
// flowx(0, y) = 5.f * std::sin(10.f * time) / N; // flowx(0, y) = 5.f * std::sin(10.f * time) / N;
// flowx(0, y) = 1.f / N; // flowx(0, y) = 10.f / N;
// flowx(N, y) = 1.f / N; // flowx(N, y) = 1.f / N;
// flowx(0, y) = 10.f * ((y + 0.5f) / N) / N; // flowx(0, y) = 10.f * ((y + 0.5f) / N) / N;
@ -175,12 +198,12 @@ void water_2d_app::update()
// if (y < N / 4) flowx(0, y) = 3.f / N; // if (y < N / 4) flowx(0, y) = 3.f / N;
// if (y >= 3 * N / 4) flowx(N, y) = -3.f / N; // if (y >= 3 * N / 4) flowx(N, y) = -3.f / N;
if (bed(0, y) < 0.75f) // if (bed(0, y) < 0.75f)
flowx(0, y) = 4.f / N; // flowx(0, y) = 4.f / N;
// flowx(0, y) = 5.f * math::sqr(std::max(0.f, std::sin(15.f * time))) / N; // flowx(0, y) = 5.f * math::sqr(std::max(0.f, std::sin(15.f * time))) / N;
if (bed(N - 1, y) < 0.75f) // if (bed(N - 1, y) < 0.75f)
flowx(N, y) = 5.f / N; // flowx(N, y) = 5.f / N;
} }
// Rain :) // Rain :)
@ -195,9 +218,13 @@ void water_2d_app::update()
// Update X flows // Update X flows
for (int y = 0; y < N; ++y) for (int y = 0; y < N; ++y)
{
for (int x = 1; x < N; ++x) for (int x = 1; x < N; ++x)
flowx(x, y) = friction * flowx(x, y) + g * dt * (water(x - 1, y) + bed(x - 1, y) - water(x, y) - bed(x, y)); flowx(x, y) = friction * flowx(x, y) + g * dt * (water(x - 1, y) + bed(x - 1, y) - water(x, y) - bed(x, y));
flowx(0, y) = friction * flowx(0, y) + g * dt * (water(N - 1, y) + bed(N - 1, y) - water(0, y) - bed(0, y));
}
// Update Y flows // Update Y flows
for (int y = 1; y < N; ++y) for (int y = 1; y < N; ++y)
for (int x = 0; x < N; ++x) for (int x = 0; x < N; ++x)
@ -236,7 +263,7 @@ void water_2d_app::update()
{ {
float outflow = 0.f; float outflow = 0.f;
outflow += std::max(0.f, - flowx(x , y )); outflow += std::max(0.f, - flowx(x , y ));
outflow += std::max(0.f, flowx(x + 1, y )); outflow += std::max(0.f, flowx((x + 1) % N, y ));
outflow += std::max(0.f, - flowy(x , y )); outflow += std::max(0.f, - flowy(x , y ));
outflow += std::max(0.f, flowy(x , y + 1)); outflow += std::max(0.f, flowy(x , y + 1));
@ -247,7 +274,7 @@ void water_2d_app::update()
float scale = std::min(1.f, max_outflow / outflow); float scale = std::min(1.f, max_outflow / outflow);
if (flowx(x, y) < 0.f) flowx(x, y) *= scale; if (flowx(x, y) < 0.f) flowx(x, y) *= scale;
if (flowx(x + 1, y) > 0.f) flowx(x + 1, y) *= scale; if (flowx((x + 1) % N, y) > 0.f) flowx((x + 1) % N, y) *= scale;
if (flowy(x, y) < 0.f) flowy(x, y) *= scale; if (flowy(x, y) < 0.f) flowy(x, y) *= scale;
if (flowy(x, y + 1) > 0.f) flowy(x, y + 1) *= scale; if (flowy(x, y + 1) > 0.f) flowy(x, y + 1) *= scale;
@ -262,7 +289,7 @@ void water_2d_app::update()
{ {
float w_old = water(x, y); float w_old = water(x, y);
water(x, y) += dt / dx / dy * (flowx(x, y) + flowy(x, y) - flowx(x + 1, y) - flowy(x, y + 1)); water(x, y) += dt / dx / dy * (flowx(x, y) + flowy(x, y) - flowx((x + 1) % N, y) - flowy(x, y + 1));
float w_avg = (w_old + water(x, y)) / 2.f; float w_avg = (w_old + water(x, y)) / 2.f;
@ -270,7 +297,7 @@ void water_2d_app::update()
if (w_avg > 0.f) if (w_avg > 0.f)
{ {
v[0] = (flowx(x, y) + flowx(x + 1, y)) / 2.f / dx / w_avg; v[0] = (flowx(x, y) + flowx((x + 1) % N, y)) / 2.f / dx / w_avg;
v[1] = (flowy(x, y) + flowy(x, y + 1)) / 2.f / dy / w_avg; v[1] = (flowy(x, y) + flowy(x, y + 1)) / 2.f / dy / w_avg;
} }
@ -278,6 +305,20 @@ void water_2d_app::update()
} }
} }
// Add coriolis force
for (int y = 0; y < N; ++y)
{
for (int x = 0; x < N; ++x)
{
auto n = math::ort(velocity(x, y)) * coriolis * dt / 2.f * std::sin(float(math::pi) / 2.f * ((y + 0.5f) * 2.f / N - 1.f));
flowx(x, y) += n[0];
flowx((x + 1) % N, y) += n[0];
flowy(x, y) += n[1];
flowy(x, y + 1) += n[1];
}
}
if (erosion_on) if (erosion_on)
{ {
// Erosion-deposition // Erosion-deposition
@ -332,6 +373,36 @@ void water_2d_app::update()
std::swap(sediment, new_sediment); std::swap(sediment, new_sediment);
} }
// Init particles
for (int i = 0; i < particles_per_frame; ++i)
if (particles.size() < max_particles)
{
math::point pos{random::uniform<float>(rng, 0.f, N), random::uniform<float>(rng, 0.f, N)};
if (water(std::floor(pos[0]), std::floor(pos[1])) > 1e-3f)
particles.push_back({pos, 0});
}
// Update particles
for (auto & p : particles)
{
p.position += velocity(std::min<int>(N - 1, std::floor(p.position[0])), std::min<int>(N - 1, std::floor(p.position[1]))) * (N * dt);
p.position[0] = math::fmod(p.position[0], float(N));
p.position[1] = math::clamp(p.position[1], {0.f, N});
p.lifetime += 1;
}
// Destroy particles
for (int i = 0; i < particles.size();)
{
if (particles[i].lifetime >= max_particle_lifetime)
{
std::swap(particles[i], particles.back());
particles.pop_back();
}
else
++i;
}
time += dt; time += dt;
} }
@ -371,8 +442,8 @@ void water_2d_app::present()
float w = 1.f - std::exp(- 1.f * water(x, y)); float w = 1.f - std::exp(- 1.f * water(x, y));
float s = 1.f - std::exp(- 1.f * sediment(x, y)); float s = 1.f - std::exp(- 1.f * sediment(x, y));
if (water(x, y) > 1e-3f) // if (water(x, y) > 1e-3f)
w = 0.5f + 0.5f * w; // w = 0.5f + 0.5f * w;
auto bed_color = gfx::to_coloru8(gfx::color_4f{0.96f, 0.72f, 0.53f, b}); auto bed_color = gfx::to_coloru8(gfx::color_4f{0.96f, 0.72f, 0.53f, b});
// auto bed_color = gfx::to_coloru8(gfx::color_4f{0.3f, 0.5f, 0.1f, b}); // auto bed_color = gfx::to_coloru8(gfx::color_4f{0.3f, 0.5f, 0.1f, b});
@ -398,20 +469,33 @@ void water_2d_app::present()
{ {
auto center = simulation_area.corner((x + 0.5f) * invN, (y + 0.5f) * invN); auto center = simulation_area.corner((x + 0.5f) * invN, (y + 0.5f) * invN);
auto v = water(x, y) * velocity(x, y) / 100.f; auto v = water(x, y) * velocity(x, y) / 40.f;
float max_length = 0.02f; float max_length = 0.05f;
auto l = math::length(v); auto l = math::length(v);
auto color = gfx::to_coloru8(gfx::color_4f{1.f, std::exp(- 10.f * l), 0.f, 1.f}); // auto color = gfx::to_coloru8(gfx::color_4f{1.f, std::exp(- 10.f * l), 0.f, 1.f});
auto color = gfx::color_rgba{255, 255, 255, 255};
float s = std::min(1.f, l / max_length); float s = std::min(1.f, l / max_length);
// s = 1.f;
painter.line(center, center + math::normalized(v) * max_length * s, s * 0.004f, s * 0.002f, color, color, false); auto d = math::normalized(v) * max_length * s / 2.f;
painter.line(center - d, center + d, s * 0.001f, 0.f, color, color, false);
} }
} }
} }
if (show_particles)
{
for (auto const & p : particles)
{
auto pos = math::lerp(simulation_area, math::vector{p.position[0] / N, p.position[1] / N});
painter.circle(pos, 0.005f, {255, 255, 255, 127}, 6);
}
}
auto mouse = math::cast<float>(state().mouse); auto mouse = math::cast<float>(state().mouse);
mouse[0] = math::lerp(view_area[0], mouse[0] / state().size[0]); mouse[0] = math::lerp(view_area[0], mouse[0] / state().size[0]);
@ -447,6 +531,11 @@ void water_2d_app::present()
} }
painter.render(transform); painter.render(transform);
{
painter.text({20.f, 20.f}, std::format("{} particles", particles.size()), {.scale = {2.f, 2.f}, .x = gfx::painter::x_align::left, .y = gfx::painter::y_align::top, .c = {0, 0, 0, 255}});
painter.render(math::window_camera{state().size[0], state().size[1]}.transform());
}
} }
void water_2d_app::on_event(app::resize_event const & event) void water_2d_app::on_event(app::resize_event const & event)
@ -467,6 +556,9 @@ void water_2d_app::on_event(app::key_event const & event)
if (event.down && event.key == app::keycode::V) if (event.down && event.key == app::keycode::V)
show_velocity ^= true; show_velocity ^= true;
if (event.down && event.key == app::keycode::P)
show_particles ^= true;
if (event.down && event.key == app::keycode::E) if (event.down && event.key == app::keycode::E)
erosion_on ^= true; erosion_on ^= true;