fixed simulation speed consistency issues, added SCRAM, added text to

buttons

assets/shader/main.vsh

@@ -8,6 +8,7 @@ layout (location = 2) in vec4 aPos;
 layout (location = 3) in vec3 aNormal;
 
 uniform mat4 model;
+uniform mat4 camera;
 uniform mat4 projection;
 
 out float brightness;
@@ -16,8 +17,8 @@ out vec2 texPos;
 
 void main()
 {
-	vec4 pos = model * aPos;
-	vec3 cNormal = vec3(0.f, 0.f, 1.f) * mat3(model);
+	vec4 pos = camera * model * aPos;
+	vec3 cNormal = vec3(0.f, 0.f, 1.f) * mat3(camera * model);
 
 	brightness = dot(normalize(aNormal), normalize(cNormal)) * 0.25f + 0.75f;
 

src/graphics/input/keyboard.cpp

@@ -9,6 +9,7 @@
 #include "../window.hpp"
 #include "../resize.hpp"
 #include "../camera.hpp"
+#include "../../system.hpp"
 
 using namespace sim::graphics;
 
@@ -24,6 +25,25 @@ static void cb_keypress(GLFWwindow* win, int key, int sc, int action, int mods)
 	if(action == GLFW_PRESS)
 	{
 		pressed[key] = true;
+
+		switch(key)
+		{
+		case GLFW_KEY_1:
+			sim::system::active.speed = 1;
+			break;
+		case GLFW_KEY_2:
+			sim::system::active.speed = 10;
+			break;
+		case GLFW_KEY_3:
+			sim::system::active.speed = 60;
+			break;
+		case GLFW_KEY_4:
+			sim::system::active.speed = 600;
+			break;
+		case GLFW_KEY_5:
+			sim::system::active.speed = 3600;
+			break;
+		}
 	}
 
 	if(action == GLFW_RELEASE)

src/graphics/mesh/glmesh.cpp

@@ -59,8 +59,7 @@ void glmesh::bind()
 
 void glmesh::uniform()
 {
-	glm::mat4 m = camera::get_matrix() * model_matrix;
-	glUniformMatrix4fv(shader::gl_model, 1, false, &m[0][0]);
+	glUniformMatrix4fv(shader::gl_model, 1, false, &model_matrix[0][0]);
 	glUniformMatrix4fv(shader::gl_tex_mat, 1, false, &colour_matrix[0][0]);
 }
 

src/graphics/mesh/mesh.cpp

@@ -101,6 +101,11 @@ bool mesh::check_focus(double len) const
 	return focus::is_triggered() && check_intersect(camera::get_pos(), camera::get_normal() * len);
 }
 
+bool mesh::check_focus() const
+{
+	return check_focus(2.5);
+}
+
 bool mesh::check_intersect(vec3 pos, vec3 path) const
 {
 	double l = glm::length(path);

src/graphics/mesh/mesh.hpp

@@ -28,6 +28,7 @@ struct mesh
 	void add(const mesh& o, glm::mat4 mat);
 
 	mesh to_lines() const;
+	bool check_focus() const;
 	bool check_focus(double len) const;
 	bool check_intersect(glm::vec<3, double> pos, glm::vec<3, double> path) const;
 	glm::vec<3, double> calc_intersect(glm::vec<3, double> pos, glm::vec<3, double> path) const;

src/graphics/monitor/core.cpp

@@ -19,7 +19,7 @@ struct core_focus_t : public focus::focus_t
 {
 	virtual void on_cursor_pos(double x, double y)
 	{
-		sim::system::active.reactor->rod_speed -= y * 1e-6;
+		sim::system::active.reactor->add_rod_speed(-y * 1e-6);
 	}
 
 	void set_all(bool state)
@@ -28,9 +28,9 @@ struct core_focus_t : public focus::focus_t
 		{
 			sim::reactor::rod* r = sim::system::active.reactor->rods[i].get();
 
-			if(r->should_select() && (r->is_selected() != state))
+			if(r->should_select())
 			{
-				r->toggle_selected();
+				r->selected = state;
 			}
 		}
 	}
@@ -77,7 +77,7 @@ struct core_focus_t : public focus::focus_t
 			toggle_auto();
 			break;
 		case GLFW_KEY_X:
-			sim::system::active.reactor->rod_speed = 0;
+			sim::system::active.reactor->reset_rod_speed();
 			break;
 		case GLFW_KEY_Q:
 			set_all(true);
@@ -112,14 +112,20 @@ void core::init()
 	mesh2.set(rmesh2, GL_STATIC_DRAW);
 	
 	mesh_click.load_model("../assets/model/", "reactor_core_input.stl");
+	mesh_scram.load_model("../assets/model/", "reactor_core_scram.stl");
 }
 
 void core::update()
 {
-	if(mesh_click.check_focus(2))
+	if(mesh_click.check_focus())
 	{
 		focus::set(std::make_unique<core_focus_t>());
 	}
+
+	if(mesh_scram.check_focus())
+	{
+		sim::system::active.reactor->scram();
+	}
 }
 
 void core::render()
@@ -169,7 +175,7 @@ void core::render()
 			mesh2.render();
 		}
 
-		if(r->is_selected())
+		if(r->selected)
 		{
 			mesh2.model_matrix = mat * mat_select;
 			mesh2.colour_matrix = arrays::colour({1, 1, 0, 1});

src/graphics/monitor/core.hpp

@@ -9,7 +9,7 @@ namespace sim::graphics::monitor
 
 class core
 {
-	sim::graphics::mesh mesh_click;
+	sim::graphics::mesh mesh_click, mesh_scram;
 	sim::graphics::glmesh mesh1, mesh2;
 	
 public:

src/graphics/monitor/vessel.cpp

@@ -37,10 +37,11 @@ void vessel::init()
 	ss << "Steam\n";
 	ss << "Pressure\n";
 	ss << "Level\n";
-	ss << "Void Ratio\n\n\n\n";
+	ss << "Void Ratio\n\n";
 	ss << "Reactor Core\n\n";
+	ss << "Energy Output\n";
+	ss << "Neutron Flux\n\n";
 	ss << "Temperature\nMin\nMax\n\n";
-	ss << "Neutron Flux\nSlow\nFast\n\n";
 	ss << "Control Rods\nMin\nMax\nSpeed\n";
 
 	rmesh.load_text(ss.str().c_str(), 0.04);
@@ -98,13 +99,13 @@ void vessel::update()
 	ss << show( sys.vessel->get_steam() ) << " g\n";
 	ss << show( sys.vessel->get_pressure() * 0.001 ) << " kPa\n";
 	ss << show( sys.vessel->get_level() ) << " / " << show( sys.vessel->get_volume() ) << " L\n";
-	ss << show( sys.vessel->get_void_ratio() * 100 ) << " %\n\n\n\n\n\n\n";
+	ss << show( sys.vessel->get_void_ratio() * 100 ) << " %\n\n\n\n";
+	ss << show( sys.reactor->get_energy_output() * 0.001 ) << " kW\n";
+	ss << show( sys.reactor->get_flux() ) << " n/cm2/s\n\n\n";
 	ss << show( temp_min ) << " C\n";
 	ss << show( temp_max ) << " C\n\n\n";
-	ss << ( sys.reactor->get_total(sim::reactor::rod::val_t::N_SLOW) * 1e12 ) << " pmol\n";
-	ss << ( sys.reactor->get_total(sim::reactor::rod::val_t::N_FAST) * 1e12 ) << " pmol\n\n\n";
-	ss << show( crod_min * 100 ) << " %\n";
-	ss << show( crod_max * 100 ) << " %\n";
+	ss << show( 100 - crod_max * 100 ) << " %\n";
+	ss << show( 100 - crod_min * 100 ) << " %\n";
 	ss << show( sys.reactor->rod_speed * 100, 1e6 ) << " %/s";
 	if(sys.reactor->rod_speed == 0) ss << " (Stopped)";
 	ss << "\n";

src/graphics/shader.cpp

@@ -15,6 +15,7 @@ static unsigned int prog_id;
 
 int shader::gl_tex_mat;
 int shader::gl_model;
+int shader::gl_camera;
 int shader::gl_projection;
 
 static int load_shader(const char* src, int type)
@@ -68,6 +69,7 @@ unsigned int shader::init_program()
 	
 	gl_tex_mat = glGetUniformLocation(prog_id, "tex_mat");
 	gl_model = glGetUniformLocation(prog_id, "model");
+	gl_camera = glGetUniformLocation(prog_id, "camera");
 	gl_projection = glGetUniformLocation(prog_id, "projection");
 
 	glUseProgram(prog_id);

src/graphics/shader.hpp

@@ -6,6 +6,7 @@ namespace sim::graphics::shader
 
 extern int gl_tex_mat;
 extern int gl_model;
+extern int gl_camera;
 extern int gl_projection;
 
 unsigned int init_program();

src/graphics/ui.cpp

@@ -0,0 +1,67 @@
+
+#include "ui.hpp"
+
+#include <GL/glew.h>
+#include <GLFW/glfw3.h>
+
+#include <glm/ext/matrix_transform.hpp>
+
+#include "mesh/mesh.hpp"
+#include "mesh/glmesh.hpp"
+#include "mesh/arrays.hpp"
+#include "mesh/font.hpp"
+#include "mesh/texture.hpp"
+#include "resize.hpp"
+#include "shader.hpp"
+
+#include "widget/clock.hpp"
+
+using namespace sim::graphics;
+
+static glmesh StaticMeshData;
+static widget::clock WidgetClock;
+
+void ui::init()
+{
+	mesh m;
+
+	unsigned int handle = texture::handle_white;
+	const unsigned int indices[] = {0, 1, 3, 0, 3, 2};
+	const arrays::vertex vertices[] = {
+		arrays::vertex(handle, {0, 0}, {-1, -1, 0, 1}, {0, 0, -1}),
+		arrays::vertex(handle, {0, 1}, {-1,  1, 0, 1}, {0, 0, -1}),
+		arrays::vertex(handle, {1, 0}, { 1, -1, 0, 1}, {0, 0, -1}),
+		arrays::vertex(handle, {1, 1}, { 1,  1, 0, 1}, {0, 0, -1}),
+	};
+
+	m.set_indices(indices, 6);
+	m.set_vertices(vertices, 4);
+
+	StaticMeshData.bind();
+	StaticMeshData.set(m, GL_STATIC_DRAW);
+	StaticMeshData.colour_matrix = glm::scale(glm::mat4(1), glm::vec3(1) * 0.75f);
+}
+
+void ui::update(double dt)
+{
+	WidgetClock.update(dt);
+}
+
+void ui::render()
+{
+	glClear(GL_DEPTH_BUFFER_BIT);
+
+	glm::vec2 wsize(resize::get_size() / 2);
+	
+	glm::mat4 mat_projection = glm::mat4(1);
+	glm::mat4 mat_camera = glm::scale(glm::mat4(1), glm::vec3(1.0f / wsize * glm::vec2(1, -1), -1));
+	glUniformMatrix4fv(shader::gl_projection, 1, false, &mat_projection[0][0]);
+	glUniformMatrix4fv(shader::gl_camera, 1, false, &mat_camera[0][0]);
+
+	StaticMeshData.bind();
+	StaticMeshData.uniform();
+	StaticMeshData.render();
+
+	WidgetClock.render();
+}
+

src/graphics/ui.hpp

@@ -0,0 +1,12 @@
+
+#pragma once
+
+namespace sim::graphics::ui
+{
+
+void init();
+void update(double dt);
+void render();
+
+};
+

src/graphics/widget/clock.cpp

@@ -0,0 +1,51 @@
+
+#include "clock.hpp"
+
+#include <GL/glew.h>
+#include <GLFW/glfw3.h>
+
+#include <glm/ext/matrix_transform.hpp>
+
+#include <iomanip>
+#include <sstream>
+#include <cmath>
+
+#include "../mesh/arrays.hpp"
+#include "../mesh/font.hpp"
+#include "../mesh/arrays.hpp"
+#include "../resize.hpp"
+#include "../../system.hpp"
+
+using namespace sim::graphics::widget;
+
+void clock::update(double dt)
+{
+	mesh m;
+	glm::vec2 wsize(resize::get_size() / 2);
+	std::stringstream ss;
+
+	int t_s = std::fmod(at, 60);
+	int t_m = std::fmod(at / 60, 60);
+	int t_h = std::fmod(at / 3600, 24);
+
+	ss << "Time: " << std::setfill('0') << std::setw(2) << t_h << ":";
+   	ss << std::setfill('0') << std::setw(2) << t_m << ":";
+	ss << std::setfill('0') << std::setw(2) << t_s << "\n";
+	ss << "Day: " << std::floor(at / (3600 * 24)) << "\n";
+	at += dt * sim::system::active.speed;
+
+	m.load_text(ss.str().c_str(), 20);
+	
+	data.bind();
+	data.model_matrix = glm::translate(glm::mat4(1), glm::vec3(-wsize + glm::vec2(2, 2), 0));
+	data.colour_matrix = arrays::colour({1, 1, 1, 1});
+	data.set(m, GL_DYNAMIC_DRAW);
+}
+
+void clock::render()
+{
+	data.bind();
+	data.uniform();
+	data.render();
+}
+

src/graphics/widget/clock.hpp

@@ -0,0 +1,19 @@
+
+#pragma once
+
+#include "../mesh/glmesh.hpp"
+
+namespace sim::graphics::widget
+{
+
+struct clock
+{
+	double at = 3600 * 12;
+	glmesh data;
+	
+	void update(double dt);
+	void render();
+};
+
+};
+

src/graphics/window.cpp

@@ -21,13 +21,14 @@
 #include "monitor/vessel.hpp"
 #include "monitor/core.hpp"
 #include "mesh/texture.hpp"
+#include "ui.hpp"
 
 using namespace sim::graphics;
 
 static GLFWwindow* win;
 static bool win_should_close = false;
 
-static glmesh MeshScene, MeshDebug;
+static glmesh MeshScene;
 static monitor::vessel MonitorVessel;
 static monitor::core MonitorCore;
 
@@ -81,8 +82,8 @@ void window::create()
 
 	glEnable(GL_DEBUG_OUTPUT);
 	glEnable(GL_DEBUG_OUTPUT_SYNCHRONOUS);
-	glEnable(GL_CULL_FACE);
 	glEnable(GL_DEPTH_TEST);
+	glEnable(GL_CULL_FACE);
 	glEnable(GL_BLEND);
 
 	glBlendFunc(GL_SRC_ALPHA, GL_ONE_MINUS_SRC_ALPHA);
@@ -93,18 +94,16 @@ void window::create()
 	resize::init();
 	texture::init();
 	font::init();
+	ui::init();
 
 	shader::init_program();
 
 	sim::system& sys = sim::system::active;
+	mesh m;
 
-	sys.scene.load_model("../assets", "scene-baked.glb");
+	m.load_model("../assets", "scene-baked.glb");
 	MeshScene.bind();
-	MeshScene.set(sys.scene, GL_STATIC_DRAW);
-
-//	sys.scene.load_model("../assets/model", "reactor_core_input.stl");
-//	MeshDebug.bind();
-//	MeshDebug.set(sys.scene, GL_STATIC_DRAW);
+	MeshScene.set(m, GL_STATIC_DRAW);
 
 	sys.scene.load_model("../assets/model", "scene_collisions.stl");
 //	MeshCollisionScene.bind();
@@ -117,22 +116,25 @@ void window::create()
 	glViewport(0, 0, 800, 600);
 }
 
-void window::loop()
+void window::update(double dt)
 {
 	glfwPollEvents();
 
 	MonitorCore.update();
 	MonitorVessel.update();
 
-	glm::mat4 mat_projection = glm::perspective(glm::radians(80.0f), resize::get_aspect(), 0.01f, 20.f);
+	ui::update(dt);
+}
+
+void window::render()
+{
+	glm::mat4 mat_camera = camera::get_matrix();
+	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]);
 
 	glClearColor(0, 0, 0, 1.0f);
 	glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT);
-
-//	MeshDebug.bind();
-//	MeshDebug.uniform();
-//	MeshDebug.render();
 	
 	MeshScene.bind();
 	MeshScene.uniform();
@@ -141,6 +143,8 @@ void window::loop()
 	MonitorCore.render();
 	MonitorVessel.render();
 
+	ui::render();
+
 	glfwSwapBuffers(win);
 }
 

src/graphics/window.hpp

@@ -10,7 +10,8 @@ namespace sim::graphics::window
 
 bool should_close();
 void create();
-void loop();
+void update(double dt);
+void render();
 void destroy();
 void close();
 

src/main.cpp

@@ -41,9 +41,11 @@ int main()
 		clock += passed;
 		
 		sim::system::active.update(dt);
+
 		graphics::camera::update(dt);
-		graphics::window::loop();
+		graphics::window::update(dt);
 		graphics::focus::update();
+		graphics::window::render();
 	}
 
 	graphics::window::destroy();

src/reactor/builder.cpp

@@ -7,7 +7,7 @@
 
 using namespace sim::reactor;
 
-sim::reactor::reactor sim::reactor::builder(const int W, const int H, const double CW, const double CH, fuel::fuel_rod fr, control::boron_rod br, coolant::vessel& v, const char** lines)
+sim::reactor::reactor sim::reactor::builder(const int W, const int H, const double CW, const double CH, fuel::fuel_rod fr, coolant::vessel& v, const char** lines)
 {
 	std::vector<std::unique_ptr<rod>> arr(W * H);
 	
@@ -23,7 +23,7 @@ sim::reactor::reactor sim::reactor::builder(const int W, const int H, const doub
 			r = std::make_unique<fuel::fuel_rod>(fr);
 			break;
 		case 'C':
-			r = std::make_unique<control::boron_rod>(br);
+			r = std::make_unique<control::boron_rod>(v);
 			break;
 		case 'G':
 			r = std::make_unique<control::graphite_rod>();

src/reactor/builder.hpp

@@ -13,7 +13,7 @@
 namespace sim::reactor
 {
 
-reactor builder(const int W, const int H, const double CW, const double CH, fuel::fuel_rod fr, control::boron_rod br, coolant::vessel& v, const char** lines);
+reactor builder(const int W, const int H, const double CW, const double CH, fuel::fuel_rod fr, coolant::vessel& v, const char** lines);
 
 };
 

src/reactor/control/boron_rod.cpp

@@ -9,9 +9,9 @@ constexpr double boron_density = 2340000; // g/m^3
 constexpr double boron_molar_mass = 10; // g/mol
 constexpr double boron_molar_density = boron_density / boron_molar_mass; // mol/m^3
 
-boron_rod::boron_rod(coolant::vessel& v, double max) : coolant::pipe(v), max(max)
+boron_rod::boron_rod(coolant::vessel& v) : coolant::pipe(v)
 {
-	
+
 }
 
 void boron_rod::display(std::ostream& o) const
@@ -39,15 +39,12 @@ void boron_rod::update(double secs)
 
 	update_rod(secs);
 
-	double k = (1 - absorbed) * inserted * max;
-	double m = 1 - std::pow(0.5, secs * -std::log2(1 - k));
+	double m = (1 - absorbed) * inserted;
 	double r_fast = vals[val_t::N_FAST] * m;
-	double r_slow = vals[val_t::N_SLOW] * m;
 
 	vals[val_t::N_FAST] -= r_fast;
-	vals[val_t::N_SLOW] -= r_slow;
-	absorbed += (r_fast + r_slow) / limit;
-	
+	absorbed += r_fast / limit;
+
 	update_pipe(secs);
 }
 

src/reactor/control/boron_rod.hpp

@@ -6,10 +6,8 @@
 namespace sim::reactor::control
 {
 
-class boron_rod : public sim::reactor::coolant::pipe
+class boron_rod : public coolant::pipe
 {
-	const double max;
-
 	double inserted = 1;
 	double absorbed = 0;
 	
@@ -19,8 +17,8 @@ class boron_rod : public sim::reactor::coolant::pipe
 	virtual int get_id() const { return 5; }
 
 public:
-
-	boron_rod(coolant::vessel& v, double max);
+	
+	boron_rod(coolant::vessel& v);
 
 	virtual void update(double secs);
 	void set_reactivity(double a);

src/reactor/coolant/pipe.cpp

@@ -12,7 +12,7 @@ pipe::pipe(coolant::vessel& v)
 
 double pipe::get_k(val_t type) const
 {
-	return vessel->get_level() / vessel->get_volume() * (1 - vessel->get_void_ratio()) * 0.5;
+	return vessel->get_level() / vessel->get_volume() * 0.5;
 }
 
 void pipe::update(double secs)
@@ -24,10 +24,10 @@ void pipe::update(double secs)
 void pipe::update_pipe(double secs)
 {
 	sim::reactor::reactor* r = (sim::reactor::reactor*)reactor;
-	double m = r->cell_width * r->cell_width * r->cell_height * 1e6;
+	double m_heat = r->cell_width * r->cell_width * r->cell_height * 1e6;
 
-	vals[val_t::HEAT] = vessel->add_heat(m, vals[val_t::HEAT]);
-	vals[val_t::N_SLOW] += vals[val_t::N_FAST];
+	vals[val_t::HEAT] = vessel->add_heat(m_heat, vals[val_t::HEAT]);
+	vals[val_t::N_SLOW] += vals[val_t::N_FAST] * (1 - vessel->get_void_ratio());
 	vals[val_t::N_FAST] = 0;
 }
 

src/reactor/fuel/fuel_rod.cpp

@@ -21,11 +21,18 @@ void fuel_rod::display(std::ostream& o) const
 	
 	o << "Fuel: " << (s.get_fuel() * mol) << " / " << (s.get_mass() * mol) << " mol\n";
 	o << "Efficiency: " << (s.get_efficiency() * 100) << " %\n";
-	o << "Output: " << (s.get_energy() * mol * energy_density) << " W\n";
+	o << "Output: " << get_energy_output() << " W\n";
 	o << "Iodine: " << (s.get_i_135() * mol) << " mol\n";
 	o << "Xenon: " << (s.get_xe_135() * mol) << " mol\n";
 }
 
+double fuel_rod::get_energy_output() const
+{
+	double mol = fuel_molar_density * get_volume();
+
+	return s.get_energy() * mol * energy_density;
+}
+
 static float map(float v, float imin, float imax, float omin, float omax)
 {
 	return (v - imin) * (omax - omin) / (imax - imin) + omin;
@@ -77,11 +84,12 @@ void fuel_rod::update(double secs)
 	
 	s.clear_energy();
 	s.clear_fast_neutrons();
+	s.clear_slow_neutrons();
 	s.add_slow_neutrons(vals[val_t::N_SLOW] / mol);
 	s.update(secs);
 
 	vals[val_t::HEAT] += s.get_energy() * mol * energy_density * secs;
 	vals[val_t::N_FAST] += s.get_fast_neutrons() * mol;
-	vals[val_t::N_SLOW] = 0;
+	vals[val_t::N_SLOW] = s.get_slow_neutrons() * mol;
 }
 

src/reactor/fuel/fuel_rod.hpp

@@ -17,6 +17,7 @@ class fuel_rod : public sim::reactor::rod
 	virtual bool has_sensors(val_t t) const { return true; }
 	virtual const char* get_name() const { return "Fuel"; }
 	virtual int get_id() const { return 1; }
+	virtual double get_energy_output() const;
 	virtual glm::vec4 get_colour() const;
 	
 public:

src/reactor/fuel/half_life.hpp

@@ -10,11 +10,15 @@ const double Te_135 = 19;
 const double I_135 = 23652;
 const double Xe_135 = 32904;
 const double Cs_135 = 41971608000000;
+const double U_235 = 2.22e16;
+const double U_238 = 1.41e17;
 
-constexpr double get(double secs, double hl) noexcept
+template <typename T>
+constexpr T get(T secs, T hl) noexcept
 {
 	return std::pow(0.5, secs / hl);
 }
 
+
 };
 

src/reactor/fuel/sample.cpp

@@ -6,26 +6,26 @@
 
 using namespace sim::reactor::fuel;
 
-static const double NEUTRON_BG = 1e-25;
+constexpr double NEUTRON_BG = 1e-30;
 
 sample::sample(double fuel)
 {
 	this->fuel = fuel;
+	this->u_238 = 1 - fuel;
 	this->mass = 1;
 }
 
 void sample::update(double secs)
 {
 	double m;
-
+	
 	// decay waste and extract products
 	waste.update(secs);
 	fast_neutrons += waste.extract_neutrons();
-	energy += waste.extract_energy() * (1.0 / 30.0);
+	energy += waste.extract_energy() * (1.0 / 30.0) / secs;
 
 	// decay Xe-135
-	m = half_life::get(secs, half_life::Xe_135);
-	xe_135 *= m;
+	xe_135 *= half_life::get(secs, half_life::Xe_135);
 
 	// decay I-135 into Xe-135
 	m = half_life::get(secs, half_life::I_135);

src/reactor/fuel/sample.hpp

@@ -19,6 +19,7 @@ class sample
 	double i_135 = 0;
 	double xe_135 = 0;
 	double te_135 = 0;
+	double u_238 = 0;
 	double mass = 0;
 
 	double energy = 0; // W
@@ -39,6 +40,7 @@ public:
 	constexpr double get_mass() const { return mass; }
 	constexpr double get_energy() const { return energy; }
 	constexpr double get_fast_neutrons() const { return fast_neutrons; }
+	constexpr double get_slow_neutrons() const { return slow_neutrons; }
 	constexpr double get_volume() const { return mass + xe_135 * Xe_135_M; }
 	constexpr double get_efficiency() const { return efficiency; }
 	constexpr double get_te_135() const { return te_135; }
@@ -47,6 +49,7 @@ public:
 	
 	constexpr void clear_energy() { energy = 0; }
 	constexpr void clear_fast_neutrons() { fast_neutrons = 0; }
+	constexpr void clear_slow_neutrons() { slow_neutrons = 0; }
 	constexpr void add_slow_neutrons(double a) { slow_neutrons += a; }
 
 	friend std::ostream& operator<<(std::ostream& o, const sample& s)

src/reactor/reactor.cpp

@@ -41,9 +41,45 @@ reactor::reactor(const reactor& o) : cell_width(o.cell_width), cell_height(o.cel
 	}
 }
 
+void reactor::reset_rod_speed()
+{
+	rod_speed = 0;
+}
+
+void reactor::add_rod_speed(double a)
+{
+	rod_speed -= a;
+
+	if(rod_speed < -0.2)
+		rod_speed = -0.2;
+	if(rod_speed > 0.2)
+		rod_speed = 0.2;
+}
+
+void reactor::scram()
+{
+	rod_speed = 0.2;
+
+	for(int i = 0; i < size; i++)
+	{
+		if(rods[i]->should_select())
+		{
+			rods[i]->selected = true;
+		}
+	}
+}
+
 void reactor::update(double secs)
 {
 	int rods_lookup[size];
+	double temp_min, temp_max;
+
+	get_stats(rod::val_t::HEAT, temp_min, temp_max);
+
+	if(temp_max > 360)
+	{
+		scram();
+	}
 
 	for(int i = 0; i < size; i++)
 	{
@@ -74,7 +110,7 @@ void reactor::update_selected(double dt)
 	{
 		rod* r = rods[i].get();
 
-		if(r->is_selected())
+		if(r->selected)
 		{
 			r->update_selected(rod_speed * dt);
 		}
@@ -155,6 +191,35 @@ double reactor::get_total(rod::val_t type)
 	return v;
 }
 
+double reactor::get_average(rod::val_t type)
+{
+	return get_total(type) / size;
+}
+
+double reactor::get_flux()
+{
+	double v = 0;
+
+	for(int i = 0; i < size; i++)
+	{
+		v += rods[i]->get_flux();
+	}
+	
+	return v / size;
+}
+
+double reactor::get_energy_output()
+{
+	double v = 0;
+
+	for(int i = 0; i < size; i++)
+	{
+		v += rods[i]->get_energy_output();
+	}
+
+	return v;
+}
+
 void reactor::get_stats(rod::val_t type, double& min, double& max)
 {
 	min = INFINITY;

src/reactor/reactor.hpp

@@ -26,9 +26,15 @@ struct reactor
 	reactor(const reactor& r);
 	reactor(reactor&& r);
 	
+	void scram();
+	void reset_rod_speed();
+	void add_rod_speed(double a);
 	void update(double secs);
 	void get_stats(rod::val_t type, double& min, double& max);
 	void get_rod_stats(int type, double& min, double& max);
+	double get_flux();
+	double get_energy_output();
+	double get_average(rod::val_t type);
 	double get_total(rod::val_t type);
 	int move_cursor(int d);
 	void toggle_selected();

src/reactor/rod.cpp

@@ -6,6 +6,9 @@
 
 using namespace sim::reactor;
 
+// Avogadro's Number
+static double N_a = 6.02214076e23;
+
 double rod::get(val_t type) const
 {
 	return vals[type];
@@ -29,8 +32,8 @@ double rod::extract(val_t type, double s, double k, double o)
 	}
 
 	double v = m * 0.5 * (get(type) - o);
-
 	vals[type] -= v;
+
 	return v;
 }
 
@@ -38,22 +41,44 @@ void rod::interact(rod* o, double secs)
 {
 	for(int i = 0; i < rod::VAL_N; i++)
 	{
-		val_t v = (val_t)i;
-		add(v, o->extract(v, secs, get_k(v), get(v)));
+		val_t t = (val_t)i;
+		double v = o->extract(t, secs, get_k(t), get(t));
+		add(t, v);
+		
+		double v2 = std::abs(v / secs);
+		o->vals_n[t] += v2;
+		vals_n[t] += v2;
 	}
 }
 
+double rod::get_flux() const
+{
+	return (vals_n[val_t::N_FAST] + vals_n[val_t::N_SLOW]) * N_a / (get_side_area() * 10000) / 4;
+}
+
 double rod::get_volume() const
 {
 	auto r = (sim::reactor::reactor*)reactor;
 	return r->cell_width * r->cell_width * r->cell_height;
 }
 
+double rod::get_side_area() const
+{
+	auto r = (sim::reactor::reactor*)reactor;
+	return r->cell_width * r->cell_height;
+}
+
 void rod::update_rod(double secs)
 {
 	// decay the free neutrons
 	double m = std::pow(0.5, secs / 879.4);
 	vals[val_t::N_FAST] *= m;
 	vals[val_t::N_SLOW] *= m;
+
+	// clear data
+	for(int i = 0; i < rod::VAL_N; i++)
+	{
+		vals_n[(val_t)i] = 0;
+	}
 }
 

src/reactor/rod.hpp

@@ -12,14 +12,15 @@ class rod
 {
 public:
 	
+	bool selected = false;
 	void* reactor = nullptr;
-	static const int VAL_N = 3;
+	static const int VAL_N = 4;
 	
 	enum val_t
 	{
 		HEAT = 0,
 		N_SLOW = 1,
-		N_FAST = 2
+		N_FAST = 2,
 	};
 
 	virtual ~rod() {};
@@ -30,16 +31,19 @@ public:
 	virtual double get(val_t type) const;
 	virtual std::unique_ptr<rod> clone() const { return std::make_unique<rod>(*this); }
 	virtual glm::vec4 get_colour() const { return {0, 0, 0, 0}; }
+	virtual double get_energy_output() const { return 0; }
 	virtual int get_id() const { return 0; }
 
 	virtual bool has_sensors(val_t t) const { return false; }
 	virtual bool should_display() const { return false; }
 	virtual bool should_select() const { return false; }
 	virtual void update_selected(double a) { }
+
+	double get_flux() const;
+	double get_side_area() const;
 	double get_volume() const;
 	
 	constexpr void toggle_selected() { selected = !selected; }
-	constexpr bool is_selected() const { return selected; }
 
 	friend std::ostream& operator<<(std::ostream& o, const rod& r)
 	{
@@ -59,7 +63,7 @@ public:
 protected:
 
 	double vals[VAL_N] = {0};
-	bool selected = false;
+	double vals_n[VAL_N] = {0};
 
 	virtual void display(std::ostream& o) const { };
 	virtual double get_k(val_t type) const { return 0; }

src/system.cpp

@@ -36,10 +36,7 @@ system::system()
 	};
 	
 	vessel = std::make_unique<reactor::coolant::vessel>(8, 10, 300, sim::coolant::WATER);
-	reactor = std::make_unique<sim::reactor::reactor>(sim::reactor::builder(19, 19, 0.25, 8,
-		reactor::fuel::fuel_rod(0.5),
-		reactor::control::boron_rod(*vessel.get(), 1),
-		*vessel.get(), layout));
+	reactor = std::make_unique<sim::reactor::reactor>(sim::reactor::builder(19, 19, 1.0 / 4.0, 4, reactor::fuel::fuel_rod(0.5), *vessel.get(), layout));
 	
 	valve = std::make_unique<coolant::valve<reactor::coolant::vessel>>(*vessel.get(), 1, 500);
 	pump = std::make_unique<coolant::pump<reactor::coolant::vessel>>(*vessel.get(), 1e4, 15);
@@ -55,7 +52,7 @@ system::system(system&& o)
 
 void system::update(double dt)
 {
-	dt *= 10;
+	dt *= speed;
 	vessel->update(dt);
 	reactor->update(dt);
 	valve->update(dt);

src/system.hpp

@@ -21,6 +21,7 @@ struct system
 	std::unique_ptr<sim::coolant::valve<sim::reactor::coolant::vessel>> valve;
 	std::unique_ptr<sim::coolant::pump<sim::reactor::coolant::vessel>> pump;
 	sim::graphics::mesh scene;
+	double speed = 1;
 
 	system();
 	system(system&& o);