Gravity example experiments
This commit is contained in:
parent
72bbdec7d9
commit
ac14483916
1 changed files with 416 additions and 33 deletions
|
|
@ -6,6 +6,11 @@
|
||||||
#include <psemek/geom/camera.hpp>
|
#include <psemek/geom/camera.hpp>
|
||||||
#include <psemek/geom/constants.hpp>
|
#include <psemek/geom/constants.hpp>
|
||||||
#include <psemek/util/clock.hpp>
|
#include <psemek/util/clock.hpp>
|
||||||
|
#include <psemek/audio/engine.hpp>
|
||||||
|
#include <psemek/audio/constants.hpp>
|
||||||
|
#include <psemek/audio/oscillator.hpp>
|
||||||
|
#include <psemek/audio/effect/compressor.hpp>
|
||||||
|
#include <psemek/prof/profiler.hpp>
|
||||||
|
|
||||||
#include <random>
|
#include <random>
|
||||||
|
|
||||||
|
|
@ -47,44 +52,197 @@ struct particle
|
||||||
geom::point<float, 2> pos;
|
geom::point<float, 2> pos;
|
||||||
geom::vector<float, 2> vel;
|
geom::vector<float, 2> vel;
|
||||||
|
|
||||||
|
float angle;
|
||||||
|
float angle_vel;
|
||||||
|
|
||||||
float radius;
|
float radius;
|
||||||
float mass;
|
float mass;
|
||||||
|
float density;
|
||||||
|
|
||||||
|
geom::vector<float, 2> delta_pos{0.f, 0.f};
|
||||||
|
geom::vector<float, 2> delta_vel{0.f, 0.f};
|
||||||
|
|
||||||
|
geom::vector<float, 2> acc{0.f, 0.f};
|
||||||
|
|
||||||
float T = 0.f;
|
float T = 0.f;
|
||||||
};
|
};
|
||||||
|
|
||||||
|
float const G = 50.f;
|
||||||
|
float const GG = 0*1000.f;
|
||||||
|
float const GC = 0*100000.f;
|
||||||
|
|
||||||
|
float const dt = 0.01f;
|
||||||
|
|
||||||
|
float const world_size = 10000.f;
|
||||||
|
geom::point world_center{0.f, 0.f};
|
||||||
|
|
||||||
|
struct sound_stream
|
||||||
|
: audio::stream
|
||||||
|
{
|
||||||
|
sound_stream(std::vector<particle> const & particles)
|
||||||
|
: particles_(particles)
|
||||||
|
{}
|
||||||
|
|
||||||
|
std::optional<std::size_t> length() const override
|
||||||
|
{
|
||||||
|
return std::nullopt;
|
||||||
|
}
|
||||||
|
|
||||||
|
std::size_t read(float * data, std::size_t sample_count) override
|
||||||
|
{
|
||||||
|
std::size_t time = played_.load() / 2;
|
||||||
|
|
||||||
|
if (velocity_.size() < particles_.size())
|
||||||
|
{
|
||||||
|
velocity_.resize(particles_.size(), 0.f);
|
||||||
|
amplitude_.resize(particles_.size(), 0.f);
|
||||||
|
oscillator_.resize(particles_.size(), 0.f);
|
||||||
|
}
|
||||||
|
|
||||||
|
float mean_v = 0.f;
|
||||||
|
for (std::size_t i = 0; i < particles_.size(); ++i)
|
||||||
|
{
|
||||||
|
float v = geom::length(particles_[i].vel);
|
||||||
|
velocity_[i] = geom::lerp(velocity_[i], v, 0.1f);
|
||||||
|
mean_v += velocity_[i];
|
||||||
|
}
|
||||||
|
mean_v /= particles_.size();
|
||||||
|
|
||||||
|
for (std::size_t i = 0; i < particles_.size(); ++i)
|
||||||
|
{
|
||||||
|
// oscillator_[i].frequency(500.f * velocity_[i] / mean_v);
|
||||||
|
// oscillator_[i].frequency(geom::lerp(100.f, 10000.f, velocity_[i] / 1000.f));
|
||||||
|
oscillator_[i].frequency(10.f * velocity_[i]);
|
||||||
|
// oscillator_[i].frequency(0.1f * geom::sqr(velocity_[i]));
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
for (std::size_t i = 0; i < sample_count; i += 2)
|
||||||
|
{
|
||||||
|
float v = 0.f;
|
||||||
|
for (std::size_t j = 0; j < particles_.size(); ++j)
|
||||||
|
v += oscillator_[j].next().imag();
|
||||||
|
|
||||||
|
if (i >= 2)
|
||||||
|
{
|
||||||
|
v = (v + data[i - 1]) / 2.f;
|
||||||
|
}
|
||||||
|
|
||||||
|
data[i + 0] = v;
|
||||||
|
data[i + 1] = v;
|
||||||
|
|
||||||
|
time += 1;
|
||||||
|
}
|
||||||
|
|
||||||
|
played_.fetch_add(sample_count);
|
||||||
|
|
||||||
|
return sample_count;
|
||||||
|
}
|
||||||
|
|
||||||
|
// The number of samples already played from this stream
|
||||||
|
std::size_t played() const override
|
||||||
|
{
|
||||||
|
return played_.load();
|
||||||
|
}
|
||||||
|
|
||||||
|
private:
|
||||||
|
std::vector<particle> const & particles_;
|
||||||
|
std::vector<float> velocity_;
|
||||||
|
std::vector<float> amplitude_;
|
||||||
|
std::vector<audio::oscillator> oscillator_;
|
||||||
|
std::atomic<std::size_t> played_;
|
||||||
|
};
|
||||||
|
|
||||||
struct myapp : app::app
|
struct myapp : app::app
|
||||||
{
|
{
|
||||||
myapp()
|
myapp()
|
||||||
: app("Test app")
|
: app("Test app", 4)
|
||||||
, rng_{std::random_device{}()}
|
, rng_{std::random_device{}()}
|
||||||
{
|
{
|
||||||
gl::ClearColor(1.f, 1.f, 1.f, 1.f);
|
gl::ClearColor(1.f, 1.f, 1.f, 1.f);
|
||||||
|
|
||||||
|
vsync(true);
|
||||||
|
|
||||||
std::uniform_real_distribution<float> d{-200.f, 200.f};
|
audio_.output()->stream(audio::compressor(std::make_shared<sound_stream>(particles_), audio::from_db(-4.f), 0.95f));
|
||||||
|
|
||||||
|
std::uniform_real_distribution<float> d{-50.f, 50.f};
|
||||||
std::uniform_real_distribution<float> rr{0.5f, 2.f};
|
std::uniform_real_distribution<float> rr{0.5f, 2.f};
|
||||||
|
std::uniform_real_distribution<float> rden{0.25f, 1.f};
|
||||||
|
std::uniform_real_distribution<float> ra{0.f, 2.f * geom::pi};
|
||||||
|
|
||||||
for (int i = 0; i < 1000; ++i)
|
float min_R = 0.f;
|
||||||
|
float max_R = 50.f;
|
||||||
|
std::uniform_real_distribution<float> rR{0.f, 1.f};
|
||||||
|
|
||||||
|
bool star = false;
|
||||||
|
|
||||||
|
if (star)
|
||||||
|
particles_.push_back({{0.f, 0.f}, {0.f, 0.f}, 0.f, 0.f, 10.f, geom::pi * 10000.f, 1.f});
|
||||||
|
|
||||||
|
// particles_.push_back({{-1.f, 0.f}, {0.f, -1.f}, 0.f, 0.f, 1.f, 1.f, 1.f});
|
||||||
|
// particles_.push_back({{ 1.f, 0.f}, {0.f, 1.f}, 0.f, 0.f, 1.f, 1.f, 1.f});
|
||||||
|
|
||||||
|
// float planet_R[] = {200.f, 300.f, 400.f};
|
||||||
|
// float planet_a[] = {0.f, geom::rad(120.f), geom::rad(240.f)};
|
||||||
|
|
||||||
|
// if(false)
|
||||||
|
for (int i = 0; i < 500; ++i)
|
||||||
{
|
{
|
||||||
float r = rr(rng_);
|
geom::vector v{0.f, 0.f};
|
||||||
float m = r * r * r;
|
|
||||||
|
|
||||||
float v = 0.1f;
|
|
||||||
|
|
||||||
geom::point<float, 2> p;
|
geom::point<float, 2> p;
|
||||||
|
float m;
|
||||||
|
float r;
|
||||||
|
float den;
|
||||||
while (true)
|
while (true)
|
||||||
{
|
{
|
||||||
p = {((i % 2) ? -300.f : 300.f) + d(rng_), d(rng_)};
|
// p = {((i % 2) ? -200.f : 200.f) + d(rng_), d(rng_)};
|
||||||
|
|
||||||
|
r = rr(rng_);
|
||||||
|
den = rden(rng_);
|
||||||
|
m = geom::pi * r * r * den;
|
||||||
|
|
||||||
|
auto a = ra(rng_);
|
||||||
|
float R = std::sqrt(rR(rng_)) * (max_R - min_R) + min_R;
|
||||||
|
// if (std::uniform_int_distribution<int>(0, 1)(rng_) == 0)
|
||||||
|
// R += 200.f;
|
||||||
|
|
||||||
|
// auto pl = std::uniform_int_distribution<int>(0, std::size(planet_R) - 1)(rng_);
|
||||||
|
// a = planet_a[pl];
|
||||||
|
// R = planet_R[pl];
|
||||||
|
|
||||||
|
// R += std::uniform_real_distribution<float>(-10.f, 10.f)(rng_);
|
||||||
|
// a += geom::rad(std::uniform_real_distribution<float>(-2.f, 2.f)(rng_));
|
||||||
|
|
||||||
|
p =
|
||||||
|
{
|
||||||
|
R * std::cos(a),
|
||||||
|
R * std::sin(a),
|
||||||
|
};
|
||||||
|
p += geom::vector{((i % 2) ? -1000.f : 1000.f), 0.f};
|
||||||
|
|
||||||
|
v[0] = {(i % 2) ? 100.f : -100.f};
|
||||||
|
|
||||||
if (std::all_of(particles_.begin(), particles_.end(), [&](particle const & q){ return geom::distance(q.pos, p) > q.radius + r; }))
|
if (std::all_of(particles_.begin(), particles_.end(), [&](particle const & q){ return geom::distance(q.pos, p) > q.radius + r; }))
|
||||||
break;
|
break;
|
||||||
}
|
}
|
||||||
particles_.push_back({p, {v * d(rng_), v * d(rng_)}, r, m});
|
particles_.push_back({p, v, 0.f, 0.f, r, m, den});
|
||||||
// particles_.push_back({{ 400.f + d(rng), d(rng)}, {0.f, 0.f}, 1.f, 1.f});
|
}
|
||||||
|
|
||||||
// particles_.push_back({{d(rng), d(rng)}, {0.f, 0.f}, 1.f, 1.f});
|
float total_M = 0.f;
|
||||||
// particles_.push_back({{d(rng), d(rng)}, {0.f, 0.f}, 1.f, 1.f});
|
if (star)
|
||||||
|
total_M += particles_[0].mass;
|
||||||
|
for (std::size_t i = star ? 1 : 0; i < particles_.size(); ++i)
|
||||||
|
total_M += particles_[i].mass / 2.f;
|
||||||
|
|
||||||
|
if(false)
|
||||||
|
for (std::size_t i = star ? 1 : 0; i < particles_.size(); ++i)
|
||||||
|
{
|
||||||
|
auto r = particles_[i].pos - geom::point{0.f, 0.f};
|
||||||
|
auto R = geom::length(r);
|
||||||
|
float V = std::sqrt(G * total_M / R);
|
||||||
|
particles_[i].vel = geom::ort(r / R) * V;
|
||||||
|
(void)V;
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
|
|
@ -112,15 +270,27 @@ struct myapp : app::app
|
||||||
|
|
||||||
geom::point<int, 2> const screen_center { window_size_[0] / 2, window_size_[1] / 2 };
|
geom::point<int, 2> const screen_center { window_size_[0] / 2, window_size_[1] / 2 };
|
||||||
auto target = camera_center_ + flip_y(geom::cast<float>(*mouse_ - screen_center) * scale);
|
auto target = camera_center_ + flip_y(geom::cast<float>(*mouse_ - screen_center) * scale);
|
||||||
|
force_target_ = target;
|
||||||
|
|
||||||
geom::point<float, 2> pos{100.f, 0.f};
|
// geom::point<float, 2> pos{100.f, 0.f};
|
||||||
|
|
||||||
geom::vector<float, 2> vel = geom::normalized(target - pos) * 500.f;
|
// geom::vector<float, 2> vel = geom::normalized(target - pos) * 40.f;
|
||||||
|
|
||||||
particles_.push_back({pos, vel, 1.f, 1.f});
|
// particles_.push_back({pos, vel, 1.f, 1.f});
|
||||||
|
|
||||||
|
// for (auto & p : particles_)
|
||||||
|
// {
|
||||||
|
// auto r = p.pos - target;
|
||||||
|
// p.vel += 4000.f * r / geom::length_sqr(r) / p.mass;
|
||||||
|
// }
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
|
void on_left_button_up() override
|
||||||
|
{
|
||||||
|
force_target_ = std::nullopt;
|
||||||
|
}
|
||||||
|
|
||||||
void on_right_button_down() override
|
void on_right_button_down() override
|
||||||
{
|
{
|
||||||
camera_drag_ = mouse_;
|
camera_drag_ = mouse_;
|
||||||
|
|
@ -144,48 +314,147 @@ struct myapp : app::app
|
||||||
camera_center_ += flip_y(geom::cast<float>(*camera_drag_ - *mouse_) * scale);
|
camera_center_ += flip_y(geom::cast<float>(*camera_drag_ - *mouse_) * scale);
|
||||||
camera_drag_ = mouse_;
|
camera_drag_ = mouse_;
|
||||||
}
|
}
|
||||||
|
|
||||||
|
if (force_target_)
|
||||||
|
{
|
||||||
|
geom::scale<float, 2> const flip_y({1.f, -1.f});
|
||||||
|
float const scale = camera_size_ / window_size_[1];
|
||||||
|
geom::point<int, 2> const screen_center { window_size_[0] / 2, window_size_[1] / 2 };
|
||||||
|
auto target = camera_center_ + flip_y(geom::cast<float>(*mouse_ - screen_center) * scale);
|
||||||
|
force_target_ = target;
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
void on_key_down(SDL_Keycode key) override
|
||||||
|
{
|
||||||
|
app::app::on_key_down(key);
|
||||||
|
|
||||||
|
if (key == SDLK_SPACE)
|
||||||
|
{
|
||||||
|
particles_.push_back({{200.f, 20.f}, {-1000.f, 0.f}, 0.f, 0.f, 1.f, 100.f, 1.f});
|
||||||
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
void update() override
|
void update() override
|
||||||
{
|
{
|
||||||
// std::shuffle(particles_.begin(), particles_.end(), rng_);
|
|
||||||
|
|
||||||
float const dt = 0.01f;
|
|
||||||
float const G = 100.f;
|
|
||||||
|
|
||||||
for (std::size_t step = 0; step < 1; ++step)
|
for (std::size_t step = 0; step < 1; ++step)
|
||||||
{
|
{
|
||||||
{
|
{
|
||||||
util::clock<> clock;
|
util::clock<> clock;
|
||||||
|
|
||||||
|
for (auto & p : particles_)
|
||||||
|
p.acc = {0.f, 0.f};
|
||||||
|
|
||||||
|
for (auto & p : particles_)
|
||||||
|
p.acc += geom::vector{0.f, -GG};
|
||||||
|
|
||||||
|
for (auto & p : particles_)
|
||||||
|
{
|
||||||
|
auto r = (p.pos - p.pos.zero());
|
||||||
|
p.acc -= GC * r / std::pow(1.f + geom::length(r), 3.f);
|
||||||
|
}
|
||||||
|
|
||||||
|
// for (std::size_t i = 0; i < particles_.size(); ++i)
|
||||||
|
// {
|
||||||
|
// log::info() << "Start: #" << i << " pos = " << std::setprecision(10) << particles_[i].pos << ", vel = " << particles_[i].vel
|
||||||
|
// << ", |vel| = " << geom::length(particles_[i].vel);
|
||||||
|
// }
|
||||||
|
|
||||||
for (std::size_t i = 0; i < particles_.size(); ++i)
|
for (std::size_t i = 0; i < particles_.size(); ++i)
|
||||||
{
|
{
|
||||||
for (std::size_t j = i + 1; j < particles_.size(); ++j)
|
for (std::size_t j = i + 1; j < particles_.size(); ++j)
|
||||||
{
|
{
|
||||||
auto const r = particles_[i].pos - particles_[j].pos;
|
auto const r = particles_[i].pos - particles_[j].pos;
|
||||||
|
|
||||||
|
// float const l = std::max(particles_[i].radius + particles_[j].radius, length(r));
|
||||||
float const l = length(r);
|
float const l = length(r);
|
||||||
|
|
||||||
auto const f = G * particles_[i].mass * particles_[j].mass * r / std::pow(l, 3.f);
|
auto const f = G * particles_[i].mass * particles_[j].mass * r / std::pow(l, 3.f);
|
||||||
|
|
||||||
particles_[i].vel -= f * dt / particles_[i].mass;
|
// log::info() << "Force: #" << i << "," << j << " = " << std::setprecision(10) << f;
|
||||||
particles_[j].vel += f * dt / particles_[j].mass;
|
|
||||||
|
particles_[i].acc -= f / particles_[i].mass;
|
||||||
|
particles_[j].acc += f / particles_[j].mass;
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
total_forces_ += clock.count();
|
total_forces_ += clock.count();
|
||||||
}
|
}
|
||||||
|
|
||||||
|
// for (std::size_t i = 0; i < particles_.size(); ++i)
|
||||||
|
// {
|
||||||
|
// log::info() << "Force: #" << i << " = " << std::setprecision(10) << particles_[i].acc;
|
||||||
|
// }
|
||||||
|
|
||||||
|
if (force_target_ && false)
|
||||||
|
{
|
||||||
|
for (auto & p : particles_)
|
||||||
|
{
|
||||||
|
auto r = p.pos - *force_target_;
|
||||||
|
p.acc += 100000.f * r / geom::length_sqr(r) / p.mass;
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
// for (auto & p : particles_)
|
||||||
|
// {
|
||||||
|
// auto old = p.pos;
|
||||||
|
// p.pos += (p.pos - p.old_pos) + p.acc * dt * dt;
|
||||||
|
// p.old_pos = old;
|
||||||
|
// }
|
||||||
|
|
||||||
|
// for (std::size_t i = 0; i < particles_.size(); ++i)
|
||||||
|
// {
|
||||||
|
// log::info() << "Integrate: #" << i << " new pos = " << std::setprecision(10) << particles_[i].pos;
|
||||||
|
// }
|
||||||
|
|
||||||
|
if (false)
|
||||||
|
{
|
||||||
|
float s = 100.f;
|
||||||
|
world_center[0] += dt * std::uniform_real_distribution<float>(-s, s)(rng_);
|
||||||
|
world_center[1] += dt * std::uniform_real_distribution<float>(-s, s)(rng_);
|
||||||
|
}
|
||||||
|
|
||||||
|
if (force_target_)
|
||||||
|
world_center = *force_target_;
|
||||||
|
|
||||||
for (auto & p : particles_)
|
for (auto & p : particles_)
|
||||||
{
|
{
|
||||||
|
p.vel += p.acc * dt;
|
||||||
p.pos += p.vel * dt;
|
p.pos += p.vel * dt;
|
||||||
|
p.angle += p.angle_vel * dt;
|
||||||
}
|
}
|
||||||
|
|
||||||
{
|
{
|
||||||
util::clock<> clock;
|
util::clock<> clock;
|
||||||
|
|
||||||
for (std::size_t iteration = 0; iteration < 10; ++iteration)
|
for (std::size_t iteration = 0; iteration < 1; ++iteration)
|
||||||
{
|
{
|
||||||
|
for (auto & p : particles_)
|
||||||
|
{
|
||||||
|
auto r = p.pos - world_center;
|
||||||
|
auto l = geom::length(r);
|
||||||
|
|
||||||
|
if (l + p.radius > world_size)
|
||||||
|
{
|
||||||
|
auto n = r / l;
|
||||||
|
|
||||||
|
auto t = geom::ort(n);
|
||||||
|
|
||||||
|
auto vn = n * geom::dot(p.vel, n);
|
||||||
|
|
||||||
|
auto vt = t * geom::dot(p.vel + t * p.angle_vel * p.radius, t);
|
||||||
|
|
||||||
|
auto pt = - vt * (1.f - std::exp(- 10.f * dt));
|
||||||
|
float I = 0.5f * p.mass * geom::sqr(p.radius);
|
||||||
|
|
||||||
|
p.pos -= n * (l + p.radius - world_size);
|
||||||
|
p.vel -= 1.75f * vn;
|
||||||
|
|
||||||
|
p.vel += pt / p.mass;
|
||||||
|
p.angle_vel += geom::det(n * p.radius, pt) / I;
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
for (std::size_t i = 0; i < particles_.size(); ++i)
|
for (std::size_t i = 0; i < particles_.size(); ++i)
|
||||||
{
|
{
|
||||||
for (std::size_t j = i + 1; j < particles_.size(); ++j)
|
for (std::size_t j = i + 1; j < particles_.size(); ++j)
|
||||||
|
|
@ -225,24 +494,111 @@ struct myapp : app::app
|
||||||
|
|
||||||
auto const vij = particles_[i].vel - particles_[j].vel;
|
auto const vij = particles_[i].vel - particles_[j].vel;
|
||||||
|
|
||||||
if (l < R && dot(vij, n) < 0.f)
|
// merging collision
|
||||||
|
if (l < R && false)
|
||||||
{
|
{
|
||||||
auto const d = n * (R - l) / 2.f * 0.5f;// * (particles_[i].mass + particles_[j].mass);
|
particle p;
|
||||||
particles_[i].pos += d;
|
|
||||||
particles_[j].pos -= d;
|
|
||||||
|
|
||||||
auto dp = - n * dot(n, vij) * 2.f / (1.f / particles_[i].mass + 1.f / particles_[j].mass);
|
auto M = (particles_[i].mass + particles_[j].mass);
|
||||||
|
p.pos = particles_[i].pos + (particles_[j].pos - particles_[i].pos) * particles_[j].mass / M;
|
||||||
|
p.vel = (particles_[i].vel * particles_[i].mass + particles_[j].vel * particles_[j].mass) / M;
|
||||||
|
p.mass = M;
|
||||||
|
p.radius = std::sqrt(geom::sqr(particles_[i].radius) + geom::sqr(particles_[j].radius));
|
||||||
|
p.density = p.mass / (geom::pi * geom::sqr(p.radius));
|
||||||
|
|
||||||
particles_[i].vel += dp / particles_[i].mass;
|
particles_[i] = p;
|
||||||
particles_[j].vel -= dp / particles_[j].mass;
|
particles_.erase(particles_.begin() + j);
|
||||||
|
break;
|
||||||
|
}
|
||||||
|
|
||||||
|
// inelastic collision
|
||||||
|
// if (l < R && dot(vij, n) < 0.f)
|
||||||
|
if (l < R && false)
|
||||||
|
{
|
||||||
|
float D = R - l;
|
||||||
|
|
||||||
|
float A = 0.5f * (1.f / particles_[i].mass + 1.f / particles_[j].mass);
|
||||||
|
float B = dot(n, vij);
|
||||||
|
auto np = n * (- B / A);
|
||||||
|
// np *= 0.5f + 0.5f * std::exp(- 0.f * dt);
|
||||||
|
|
||||||
|
auto p = particles_[j].pos + (particles_[j].radius - D / 2.f) * n;
|
||||||
|
|
||||||
|
auto ri = p - particles_[i].pos;
|
||||||
|
auto rj = p - particles_[j].pos;
|
||||||
|
|
||||||
|
auto t = geom::ort(n);
|
||||||
|
|
||||||
|
auto vri = geom::dot(t, particles_[i].vel + geom::ort(ri) * particles_[i].angle_vel);
|
||||||
|
auto vrj = geom::dot(t, particles_[j].vel + geom::ort(rj) * particles_[j].angle_vel);
|
||||||
|
|
||||||
|
auto vrij = vri - vrj;
|
||||||
|
|
||||||
|
auto tp = - t * vrij * (1.f - std::exp(- 10.f * dt));;
|
||||||
|
|
||||||
|
particles_[i].vel += (np + tp) / particles_[i].mass;
|
||||||
|
particles_[j].vel -= (np + tp) / particles_[j].mass;
|
||||||
|
|
||||||
|
float Ii = 0.5f * particles_[i].mass * geom::sqr(particles_[i].radius);
|
||||||
|
float Ij = 0.5f * particles_[j].mass * geom::sqr(particles_[j].radius);
|
||||||
|
|
||||||
|
auto ti = geom::det(ri, tp);
|
||||||
|
auto tj = geom::det(rj, -tp);
|
||||||
|
|
||||||
|
particles_[i].angle_vel += ti / Ii;
|
||||||
|
particles_[j].angle_vel += tj / Ij;
|
||||||
|
|
||||||
|
float ki = D * particles_[j].mass / (particles_[i].mass + particles_[j].mass);
|
||||||
|
float kj = - D * particles_[i].mass / (particles_[i].mass + particles_[j].mass);
|
||||||
|
|
||||||
|
particles_[i].pos += ki * n;
|
||||||
|
particles_[j].pos += kj * n;
|
||||||
|
}
|
||||||
|
|
||||||
|
// collision replaced by force
|
||||||
|
if (l < R)
|
||||||
|
{
|
||||||
|
auto f = 10000.f * n * (R - l);
|
||||||
|
particles_[i].vel += dt * f / particles_[i].mass;
|
||||||
|
particles_[j].vel -= dt * f / particles_[j].mass;
|
||||||
|
|
||||||
|
auto vn = geom::dot(vij, n) * n;
|
||||||
|
// vn *= 1.f - std::exp(- (1.f - l / R));
|
||||||
|
vn *= 0.5f;
|
||||||
|
|
||||||
|
particles_[i].vel -= vn * particles_[j].mass / (particles_[i].mass + particles_[j].mass);
|
||||||
|
particles_[j].vel += vn * particles_[i].mass / (particles_[i].mass + particles_[j].mass);
|
||||||
|
}
|
||||||
|
|
||||||
|
// Verlet collision
|
||||||
|
if (l < R && false)
|
||||||
|
{
|
||||||
|
float D = (R - l) * 1.f;
|
||||||
|
float ki = D * particles_[j].mass / (particles_[i].mass + particles_[j].mass);
|
||||||
|
float kj = - D * particles_[i].mass / (particles_[i].mass + particles_[j].mass);
|
||||||
|
|
||||||
|
auto di = ki * n;
|
||||||
|
auto dj = kj * n;
|
||||||
|
|
||||||
|
particles_[i].pos += di;
|
||||||
|
particles_[j].pos += dj;
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
|
for (auto & p : particles_)
|
||||||
|
{
|
||||||
|
p.vel += p.delta_vel;
|
||||||
|
p.pos += p.delta_pos;
|
||||||
|
p.delta_pos = {0.f, 0.f};
|
||||||
|
p.delta_vel = {0.f, 0.f};
|
||||||
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
total_collisions_ += clock.count();
|
total_collisions_ += clock.count();
|
||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
float Ep = 0.f;
|
float Ep = 0.f;
|
||||||
float Ek = 0.f;
|
float Ek = 0.f;
|
||||||
|
|
||||||
|
|
@ -256,20 +612,41 @@ struct myapp : app::app
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
log::info() << "Energy: " << Ek << " - " << (-Ep) << " = " << (Ek + Ep);
|
float omega = 0.f;
|
||||||
|
for (std::size_t i = 0; i < particles_.size(); ++i)
|
||||||
|
{
|
||||||
|
omega += geom::det(particles_[i].mass * particles_[i].vel, particles_[i].pos - geom::point{0.f, 0.f});
|
||||||
|
}
|
||||||
|
|
||||||
|
// log::info() << "Angular velocity: " << omega;
|
||||||
|
|
||||||
|
// log::info() << "Energy: " << Ek << " - " << (-Ep) << " = " << (Ek + Ep);
|
||||||
|
|
||||||
++frame_count_;
|
++frame_count_;
|
||||||
|
|
||||||
|
// std::reverse(particles_.begin(), particles_.end());
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
void present() override
|
void present() override
|
||||||
{
|
{
|
||||||
|
gl::ClearColor(0.f, 0.f, 0.1f, 0.f);
|
||||||
gl::Clear(gl::COLOR_BUFFER_BIT);
|
gl::Clear(gl::COLOR_BUFFER_BIT);
|
||||||
|
|
||||||
|
painter_.circle(world_center, world_size, {255, 255, 255, 255}, 72);
|
||||||
|
|
||||||
for (auto & p : particles_)
|
for (auto & p : particles_)
|
||||||
{
|
{
|
||||||
float c = 2.f / (std::exp(-p.T / 100000.f) + 1.f) - 1.f;
|
// float c = 2.f / (std::exp(-p.T / 100000.f) + 1.f) - 1.f;
|
||||||
painter_.circle(p.pos, p.radius, {static_cast<std::uint8_t>(c * 255.f), 0, 0, 255});
|
float c = 1.f - p.density;
|
||||||
|
auto x = static_cast<std::uint8_t>(c * 255.f);
|
||||||
|
|
||||||
|
float s = window_size_[1] / camera_size_;
|
||||||
|
|
||||||
|
float r = std::max(p.radius * s, 1.5f) / s;
|
||||||
|
painter_.circle(p.pos, r, {x, x, x, 255});
|
||||||
|
painter_.line(p.pos, p.pos + r * geom::direction(p.angle), r / 16.f, {255, 0, 0, 255}, false);
|
||||||
|
// painter_.circle(p.pos, r * 0.75f, {255, 255, 255, 255});
|
||||||
}
|
}
|
||||||
|
|
||||||
geom::orthographic_camera camera;
|
geom::orthographic_camera camera;
|
||||||
|
|
@ -287,6 +664,8 @@ struct myapp : app::app
|
||||||
{
|
{
|
||||||
log::info() << "Avg forces time: " << (total_forces_ / frame_count_);
|
log::info() << "Avg forces time: " << (total_forces_ / frame_count_);
|
||||||
log::info() << "Avg collision time: " << (total_collisions_ / frame_count_);
|
log::info() << "Avg collision time: " << (total_collisions_ / frame_count_);
|
||||||
|
|
||||||
|
prof::dump();
|
||||||
}
|
}
|
||||||
|
|
||||||
private:
|
private:
|
||||||
|
|
@ -294,6 +673,8 @@ private:
|
||||||
|
|
||||||
geom::vector<int, 2> window_size_;
|
geom::vector<int, 2> window_size_;
|
||||||
|
|
||||||
|
std::optional<geom::point<float, 2>> force_target_;
|
||||||
|
|
||||||
geom::point<float, 2> camera_center_ { 0.f, 0.f };
|
geom::point<float, 2> camera_center_ { 0.f, 0.f };
|
||||||
float camera_size_ = 500.f;
|
float camera_size_ = 500.f;
|
||||||
float camera_ratio_ = 1.f;
|
float camera_ratio_ = 1.f;
|
||||||
|
|
@ -308,6 +689,8 @@ private:
|
||||||
int frame_count_ = 0;
|
int frame_count_ = 0;
|
||||||
float total_forces_ = 0.f;
|
float total_forces_ = 0.f;
|
||||||
float total_collisions_ = 0.f;
|
float total_collisions_ = 0.f;
|
||||||
|
|
||||||
|
audio::engine audio_;
|
||||||
};
|
};
|
||||||
|
|
||||||
int main()
|
int main()
|
||||||
|
|
|
||||||
Loading…
Add table
Reference in a new issue