improve grid model

This commit is contained in:
Jay Robson 2024-02-20 16:48:01 +11:00
parent d9404cec21
commit 43db2f812f
13 changed files with 331 additions and 160 deletions

View File

@ -9,11 +9,12 @@ in vec2 texPos;
out vec4 FragColour;
uniform mat4 tex_mat;
uniform vec3 brightness;
void main()
{
vec4 texdata = texture2D(tex, texPos);
FragColour = tex_mat * texdata * colour;
FragColour = tex_mat * texdata * colour * vec4(brightness, 1);
if(FragColour.a == 0) discard;
}

139
src/electric/generator.cpp Normal file
View File

@ -0,0 +1,139 @@
#include "generator.hpp"
#include "../util/math.hpp"
#include "../util/streams.hpp"
#include <cmath>
#include <iostream>
using namespace Sim::Electric;
Generator::Generator(Turbine* turbine, Grid* grid, double length, double diameter, double mass) :
length(length), diameter(diameter), mass(mass), turbine(turbine), grid(grid)
{
}
Generator::Generator(const Json::Value& node, Turbine* turbine, Grid* grid) :
length(node["length"].asDouble()),
diameter(node["diameter"].asDouble()),
friction(node["friction"].asDouble()),
mass(node["mass"].asDouble()),
turbine(turbine),
grid(grid)
{
velocity = node["velocity"].asDouble();
phase = node["phase"].asDouble();
breaker_closed = node["breaker_closed"].asBool();
energy_generated = node["energy_generated"].asDouble();
}
void Generator::update(double dt)
{
double energy_input = turbine->extract_energy();
double energy_friction = get_rpm() / 60 * dt * friction;
double work = Util::Math::j_to_ms2(energy_input - energy_friction, mass);
phase = std::fmod(phase + Util::Math::map( get_rpm(), 0, 3600, 0, 120 * M_PI ) * dt, 2 * M_PI);
// do energy transfer stuff here
if(breaker_closed)
{
double a = get_phase_diff();
double dist_extra = 0.1;
if(is_stable || (a < 1e-5 && std::abs(get_rpm() - 60 * grid->frequency) < 1e-3))
{
is_stable = true;
energy_generated = (energy_input - energy_friction) / dt;
grid->pull_energy(energy_friction - energy_input);
phase -= a;
set_rpm(grid->frequency * 60);
return;
}
glm::vec<2, double> point(std::cos(a), std::sin(a));
glm::vec<2, double> diff1 = point - glm::vec<2, double>(dist_extra + 1, 0);
glm::vec<2, double> diff2 = point - glm::vec<2, double>(-dist_extra - 1, 0);
double strength = 1e10;
double dist1 = glm::length(diff1);
double dist2 = glm::length(diff2);
double force1_len = -strength / (4 * M_PI * dist1 * dist1);
double force2_len = strength / (4 * M_PI * dist2 * dist2);
glm::vec<2, double> force1 = diff1 / dist1 * force1_len;
glm::vec<2, double> force2 = diff2 / dist2 * force2_len;
// calc the projected force vector
double t1 = std::tan(a);
double t2 = std::tan(a + M_PI/2);
glm::vec<2, double> proj1, proj2;
proj1.x = ((point.x + force1.x) * t1 - point.x * t2 - force1.y) / (t1 - t2);
proj2.x = ((point.x + force2.x) * t1 - point.x * t2 - force2.y) / (t1 - t2);
proj1.y = (proj1.x - point.x) * t2 + point.y;
proj2.y = (proj2.x - point.x) * t2 + point.y;
glm::mat<2, 2, double> rot_mat = {
point.x, -point.y,
point.y, point.x,
};
double eddy = (get_rpm() - 60 * grid->frequency) * dt * 1e5;
// calc the amount of actual work (in change in m/s) done
glm::vec<2, double> proj = rot_mat * (proj1 + proj2) * 0.5;
double work_done = proj.y / (mass * 0.001) * dt - eddy;
work += work_done;
double e = Util::Math::ms2_to_j(work_done / dt, mass);
energy_generated = -e / dt;
grid->pull_energy(e);
if(grid->is_above_limit())
{
breaker_closed = false;
}
}
else
{
energy_generated = 0;
is_stable = false;
}
velocity = std::max(velocity + work, 0.0);
}
double Generator::get_rpm() const
{
return velocity / (M_PI * mass * 0.001 * diameter * diameter * 0.25);
}
void Generator::set_rpm(double rpm)
{
velocity = rpm * (M_PI * mass * 0.001 * diameter * diameter * 0.25);
}
double Generator::get_phase_diff() const
{
return Util::Math::mod(phase - grid->get_phase() + M_PI, 2*M_PI) - M_PI;
}
Generator::operator Json::Value() const
{
Json::Value node;
node["length"] = length;
node["diameter"] = diameter;
node["velocity"] = velocity;
node["friction"] = friction;
node["breaker_closed"] = breaker_closed;
node["energy_generated"] = energy_generated;
node["phase"] = phase;
node["mass"] = mass;
return node;
}

View File

@ -0,0 +1,46 @@
#pragma once
#include "grid.hpp"
#include "turbine.hpp"
namespace Sim::Electric
{
class Generator
{
Grid* const grid;
Turbine* const turbine;
const double length; // m
const double diameter; // m
const double friction = 1e5; // J/rev
const double mass; // g
double energy_generated = 0; // W
double velocity = 0; // m/s
double phase = 0;
void set_rpm(double rpm);
public:
bool breaker_closed = false;
bool is_stable = false;
Generator(Turbine* turbine, Grid* grid, double length, double diameter, double mass);
Generator(const Json::Value& node, Turbine* turbine, Grid* grid);
void update(double dt);
double get_rpm() const;
double get_phase_diff() const;
operator Json::Value() const;
constexpr double get_energy_generated() const { return energy_generated; }
constexpr double get_phase() const { return phase; }
};
};

55
src/electric/grid.cpp Normal file
View File

@ -0,0 +1,55 @@
#include "grid.hpp"
#include <cmath>
#include <iostream>
using namespace Sim::Electric;
Grid::Grid(const Json::Value& node) :
frequency(node["frequency"].asDouble()),
power_limit(node["power_limit"].asDouble())
{
phase = node["phase"].asDouble();
energy_output = node["energy_output"].asDouble();
power_output = node["power_output"].asDouble();
}
Grid::operator Json::Value() const
{
Json::Value node;
node["phase"] = phase;
node["frequency"] = frequency;
node["energy_output"] = energy_output;
node["power_output"] = power_output;
node["power_limit"] = power_limit;
return node;
}
void Grid::update(double dt)
{
phase = std::fmod(phase + 60.0 * 2*M_PI * dt, 2*M_PI);
power_output = energy_output / dt;
energy_output = 0;
}
void Grid::pull_energy(double joules)
{
energy_output -= joules;
}
double Grid::get_light_intensity()
{
if(energy_output < 0)
{
return std::max(0.0, 1 + energy_output / 1e8);
}
return 1;
}
bool Grid::is_above_limit()
{
return std::abs(power_output) > power_limit;
}

38
src/electric/grid.hpp Normal file
View File

@ -0,0 +1,38 @@
#pragma once
#include <vector>
#include <json/json.h>
namespace Sim::Electric
{
class Grid
{
double phase = 0;
double power_output = 0; // W
double energy_output = 0; // J
public:
const double frequency = 60;
const double power_limit = 1e9; // W
const double output_resistance = 10; // Ohms
constexpr Grid() {}
Grid(const Json::Value& node);
operator Json::Value() const;
void update(double dt);
void pull_energy(double joules);
double get_light_intensity();
bool is_above_limit();
constexpr double get_phase() { return phase; }
constexpr double get_power_output() { return power_output; }
};
};

View File

@ -9,140 +9,27 @@
using namespace Sim::Electric;
constexpr static double calc_cylinder(double h, double d)
{
double r = d / 2;
return M_PI * r * r * h * 1000;
}
Turbine::Turbine(Coolant::Fluid type, Coolant::Condenser* condenser, double length, double diameter, double mass) :
length(length), diameter(diameter), condenser(condenser),
Sim::Coolant::FluidHolder(type, calc_cylinder(length, diameter), mass)
Turbine::Turbine(Coolant::Condenser* condenser) :
condenser(condenser), Sim::Coolant::FluidHolder(condenser->fluid, 1, 0)
{
}
Turbine::Turbine(const Json::Value& node, Coolant::Condenser* condenser) :
condenser(condenser),
length(node["length"].asDouble()),
diameter(node["diameter"].asDouble()),
friction(node["friction"].asDouble()),
Sim::Coolant::FluidHolder(node)
{
velocity = node["velocity"].asDouble();
phase = node["phase"].asDouble();
breaker_closed = node["breaker_closed"].asBool();
energy_input = node["energy_input"].asDouble();
energy_generated = node["energy_generated"].asDouble();
}
void Turbine::update(double dt)
{
double energy_friction = get_rpm() / 60 * dt * friction;
double work = Util::Math::j_to_ms2(energy_input - energy_friction, extra_mass);
phase = std::fmod(phase + Util::Math::map( get_rpm(), 0, 60, 0, 2 * M_PI ) * dt, 2 * M_PI);
// do energy transfer stuff here
if(breaker_closed)
{
double a = get_phase_diff();
double dist_extra = 0.1;
if(is_stable || (a < 1e-5 && std::abs(get_rpm() - 3600) < 1e-3))
{
is_stable = true;
energy_generated = (energy_input - energy_friction) / dt;
energy_input = 0;
phase -= a;
set_rpm(3600);
return;
}
glm::vec<2, double> point(std::cos(a), std::sin(a));
glm::vec<2, double> diff1 = point - glm::vec<2, double>(dist_extra + 1, 0);
glm::vec<2, double> diff2 = point - glm::vec<2, double>(-dist_extra - 1, 0);
double strength = 1e10;
double dist1 = glm::length(diff1);
double dist2 = glm::length(diff2);
double force1_len = -strength / (4 * M_PI * dist1 * dist1);
double force2_len = strength / (4 * M_PI * dist2 * dist2);
glm::vec<2, double> force1 = diff1 / dist1 * force1_len;
glm::vec<2, double> force2 = diff2 / dist2 * force2_len;
// calc the projected force vector
double t1 = std::tan(a);
double t2 = std::tan(a + M_PI/2);
glm::vec<2, double> proj1, proj2;
proj1.x = ((point.x + force1.x) * t1 - point.x * t2 - force1.y) / (t1 - t2);
proj2.x = ((point.x + force2.x) * t1 - point.x * t2 - force2.y) / (t1 - t2);
proj1.y = (proj1.x - point.x) * t2 + point.y;
proj2.y = (proj2.x - point.x) * t2 + point.y;
glm::mat<2, 2, double> rot_mat = {
point.x, -point.y,
point.y, point.x,
};
double eddy = (get_rpm() - 3600) * dt * 1e5;
// calc the amount of actual work (in change in m/s) done
glm::vec<2, double> proj = rot_mat * (proj1 + proj2) * 0.5;
double work_done = proj.y / (extra_mass * 0.001) * dt - eddy;
work += work_done;
energy_generated = -Util::Math::ms2_to_j(work_done / dt, extra_mass) / dt;
}
else
{
energy_generated = 0;
is_stable = false;
}
velocity = std::max(velocity + work, 0.0);
energy_input = 0;
}
double Turbine::get_rpm() const
{
return velocity / (M_PI * extra_mass * 0.001 * diameter * diameter * 0.25);
}
void Turbine::set_rpm(double rpm)
{
velocity = rpm * (M_PI * extra_mass * 0.001 * diameter * diameter * 0.25);
}
double Turbine::get_phase_diff() const
{
double phase_g = System::active.clock * 120 * M_PI;
return Util::Math::mod(phase - phase_g + M_PI, 2*M_PI) - M_PI;
}
void Turbine::add_gas(double steam, double air, double t)
{
energy_input += (steam + air) * fluid.jPg; // J
energy_output += (steam + air) * fluid.jPg; // J
condenser->add_gas(steam, air, t);
}
Turbine::operator Json::Value() const
void Turbine::update(double dt)
{
Json::Value node(FluidHolder::operator::Json::Value());
node["length"] = length;
node["diameter"] = diameter;
node["velocity"] = velocity;
node["friction"] = friction;
node["breaker_closed"] = breaker_closed;
node["energy_input"] = energy_input;
node["energy_generated"] = energy_generated;
node["phase"] = phase;
return node;
}
double Turbine::extract_energy()
{
double e = energy_output;
energy_output = 0;
return e;
}

View File

@ -11,27 +11,14 @@ class Turbine : public Sim::Coolant::FluidHolder
{
Coolant::Condenser* const condenser;
const double length;
const double diameter;
const double friction = 1e5; // J/rev
double energy_input = 0; // J
double energy_generated = 0; // W
double velocity = 0; // m/s
double phase = 0;
void set_rpm(double rpm);
double energy_output = 0; // J
public:
bool breaker_closed = false;
bool is_stable = false;
Turbine(Coolant::Fluid type, Coolant::Condenser* condenser, double length, double diameter, double mass);
Turbine(const Json::Value& node, Coolant::Condenser* condenser);
Turbine(Coolant::Condenser* condenser);
void update(double dt);
double get_rpm() const;
double extract_energy();
virtual double add_heat(double m, double t) { return condenser->add_heat(m, t); }
virtual double extract_fluid(double amount) { return condenser->extract_fluid(amount); }
@ -51,12 +38,6 @@ public:
virtual double get_thermal_mass() const { return condenser->get_thermal_mass(); } // grams
virtual double get_pressure() const { return condenser->get_pressure(); } // pascals
virtual double get_gas_density() const { return condenser->get_gas_density(); } // g/L
constexpr double get_energy_generated() const { return energy_generated; }
constexpr double get_phase() const { return phase; }
double get_phase_diff() const;
operator Json::Value() const;
};
};

View File

@ -69,7 +69,7 @@ void Turbine::update(double dt)
ss << "\n\n";
ss << show( sys.turbine->get_heat() ) << " C\n";
ss << show( sys.turbine->get_pressure() / 1000 ) << " kPa\n";
ss << show( sys.turbine->get_rpm() ) << " r/min\n";
ss << show( sys.generator->get_rpm() ) << " r/min\n";
rmesh2.load_text(ss.str().c_str(), 0.04);
rmesh.add(rmesh2, glm::translate(glm::mat4(1), glm::vec3(0.5, 0, 0)));
@ -77,8 +77,8 @@ void Turbine::update(double dt)
ss = std::stringstream();
ss << "Local\n\n";
ss << show( sys.turbine->get_rpm() / 60 ) << " Hz\n";
Util::Streams::show_units( ss, sys.turbine->get_energy_generated() ) << "W\n";
ss << show( sys.generator->get_rpm() / 60 ) << " Hz\n";
Util::Streams::show_units( ss, sys.generator->get_energy_generated() ) << "W\n";
rmesh2.load_text(ss.str().c_str(), 0.04);
rmesh.add(rmesh2, glm::translate(glm::mat4(1), glm::vec3(0.4, 0.7, 0)));
@ -86,7 +86,7 @@ void Turbine::update(double dt)
ss = std::stringstream();
ss << "Grid\n\n";
ss << show( 60 ) << " Hz\n";
ss << show( sys.grid->frequency ) << " Hz\n";
rmesh2.load_text(ss.str().c_str(), 0.04);
rmesh.add(rmesh2, glm::translate(glm::mat4(1), glm::vec3(0.7, 0.7, 0)));
@ -95,21 +95,21 @@ void Turbine::update(double dt)
mesh2.set(rmesh, GL_DYNAMIC_DRAW);
}
double rpm = sys.turbine->get_rpm();
double rpm = sys.generator->get_rpm();
if(rpm > 3570 && rpm < 3630)
{
glm::mat4 mat = glm::mat4(1);
mat = glm::translate(mat, glm::vec3(6.35, 3.949, 1.35));
mat = glm::rotate(mat, float(sys.turbine->get_phase_diff()), glm::vec3(0, 1, 0));
mat = glm::rotate(mat, float(sys.generator->get_phase_diff()), glm::vec3(0, 1, 0));
mat = glm::translate(mat, glm::vec3(-6.35, -3.949, -1.35));
gm_synchroscope_dial.model_matrix = mat;
}
if(m_switch_breaker.check_focus())
sys.turbine->breaker_closed = !sys.turbine->breaker_closed;
sys.generator->breaker_closed = !sys.generator->breaker_closed;
gm_switch_breaker.model_matrix = glm::translate(glm::mat4(1), glm::vec3(0, sys.turbine->breaker_closed ? 0.07 : 0, 0));
gm_switch_breaker.model_matrix = glm::translate(glm::mat4(1), glm::vec3(0, sys.generator->breaker_closed ? 0.07 : 0, 0));
}
void Turbine::render()
@ -122,7 +122,7 @@ void Turbine::render()
mesh2.uniform();
mesh2.render();
double rpm = System::active.turbine->get_rpm();
double rpm = System::active.generator->get_rpm();
if(rpm > 3570 && rpm < 3630)
{

View File

@ -17,6 +17,7 @@ int Shader::gl_tex_mat;
int Shader::gl_model;
int Shader::gl_camera;
int Shader::gl_projection;
int Shader::gl_brightness;
static int load_shader(const char* src, int type)
{
@ -71,6 +72,7 @@ unsigned int Shader::init_program()
gl_model = glGetUniformLocation(prog_id, "model");
gl_camera = glGetUniformLocation(prog_id, "camera");
gl_projection = glGetUniformLocation(prog_id, "projection");
gl_brightness = glGetUniformLocation(prog_id, "brightness");
glUseProgram(prog_id);
glDeleteShader(vsh_id);

View File

@ -8,6 +8,7 @@ extern int gl_tex_mat;
extern int gl_model;
extern int gl_camera;
extern int gl_projection;
extern int gl_brightness;
unsigned int init_program();

View File

@ -25,6 +25,7 @@
#include "monitor/secondary_loop.hpp"
#include "monitor/turbine.hpp"
#include "mesh/texture.hpp"
#include "../system.hpp"
#include "ui.hpp"
using namespace Sim::Graphics;
@ -177,9 +178,11 @@ void Window::render()
glm::mat4 mat_camera = Camera::get_matrix();
mat_camera = glm::scale(mat_camera, {1, 1, -1});
glm::vec3 brightness = glm::vec3(Sim::System::active.grid->get_light_intensity());
glm::mat4 mat_projection = glm::perspective(glm::radians(90.0f), Resize::get_aspect(), 0.01f, 20.f);
glUniformMatrix4fv(Shader::gl_projection, 1, false, &mat_projection[0][0]);
glUniformMatrix4fv(Shader::gl_camera, 1, false, &mat_camera[0][0]);
glUniform3fv(Shader::gl_brightness, 1, &brightness[0]);
projection_matrix = mat_projection;
glClearColor(0, 0, 0, 1.0f);
@ -195,6 +198,9 @@ void Window::render()
render_scene();
brightness = glm::vec3(1);
glUniform3fv(Shader::gl_brightness, 1, &brightness[0]);
UI::render();
Focus::render_ui();

View File

@ -41,7 +41,10 @@ System::System()
vessel = std::make_unique<Reactor::Coolant::Vessel>(Sim::Coolant::WATER, 8, 10, 6e6, 5e5, 10);
reactor = std::make_unique<Reactor::Reactor>(Sim::Reactor::Builder(19, 19, 1.0 / 4.0, 4, Reactor::Fuel::FuelRod(0.2), vessel.get(), layout));
condenser = std::make_unique<Coolant::Condenser>(Sim::Coolant::WATER, 6, 4, 3e6, 30000);
turbine = std::make_unique<Electric::Turbine>(Sim::Coolant::WATER, condenser.get(), 6, 3, 2e6);
grid = std::make_unique<Electric::Grid>();
turbine = std::make_unique<Electric::Turbine>(condenser.get());
generator = std::make_unique<Electric::Generator>(turbine.get(), grid.get(), 6, 3, 2e6);
sink = std::make_unique<Coolant::Sink>(Sim::Coolant::WATER, 11, 0, 0);
evaporator = std::make_unique<Coolant::Evaporator>(Sim::Coolant::WATER, 2, 30, 0, 1000);
@ -62,7 +65,10 @@ System::System(const Json::Value& node)
vessel = std::make_unique<Reactor::Coolant::Vessel>(node["vessel"]);
reactor = std::make_unique<Reactor::Reactor>(node["reactor"], vessel.get());
condenser = std::make_unique<Coolant::Condenser>(node["condenser"]);
turbine = std::make_unique<Electric::Turbine>(node["turbine"], condenser.get());
grid = std::make_unique<Electric::Grid>(node["grid"]);
turbine = std::make_unique<Electric::Turbine>(condenser.get());
generator = std::make_unique<Electric::Generator>(node["generator"], turbine.get(), grid.get());
evaporator = std::make_unique<Coolant::Evaporator>(node["evaporator"]);
sink = std::make_unique<Coolant::Sink>(evaporator->fluid, 11, 0, 0);
@ -81,13 +87,16 @@ void System::update(double dt)
dt *= speed;
clock += dt;
grid->update(dt);
reactor->update(dt);
vessel->update(dt);
turbine_inlet_valve->update(dt);
turbine_bypass_valve->update(dt);
condenser->update(dt);
turbine->update(dt);
generator->update(dt);
primary_pump->update(dt);
secondary_pump->update(dt);
freight_pump->update(dt);
@ -100,8 +109,9 @@ System::operator Json::Value() const
{
Json::Value node;
node["grid"] = *grid;
node["vessel"] = *vessel;
node["turbine"] = *turbine;
node["generator"] = *generator;
node["condenser"] = *condenser;
node["evaporator"] = *evaporator;
node["pump"]["primary"] = *primary_pump;

View File

@ -13,6 +13,8 @@
#include "coolant/evaporator.hpp"
#include "coolant/sink.hpp"
#include "electric/turbine.hpp"
#include "electric/generator.hpp"
#include "electric/grid.hpp"
namespace Sim
{
@ -28,7 +30,10 @@ struct System
std::unique_ptr<Sim::Coolant::Condenser> condenser;
std::unique_ptr<Sim::Coolant::CondenserSecondary> condenser_secondary;
std::unique_ptr<Sim::Coolant::Evaporator> evaporator;
std::unique_ptr<Sim::Electric::Turbine> turbine;
std::unique_ptr<Sim::Electric::Generator> generator;
std::unique_ptr<Sim::Electric::Grid> grid;
std::unique_ptr<Sim::Coolant::Pump> primary_pump;
std::unique_ptr<Sim::Coolant::Pump> secondary_pump;