jay
/
fast-nuclear-sim
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

added some things but its a buggy mess

This commit is contained in:
Jay Robson 2024-02-05 18:33:31 +11:00
parent 67b9f29745
commit b8191cdf65
16 changed files with 186 additions and 111 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

4
src/coolant/condenser.cpp
View File

@ -13,8 +13,8 @@ constexpr static double calc_cylinder(double h, double d)
return M_PI * r * r * h;
}
condenser::condenser(fluid_t type, double height, double diameter, double level) :
height(height), diameter(diameter), fluid_holder(type, calc_cylinder(height, diameter))
condenser::condenser(fluid_t type, double height, double diameter, double mass, double level) :
height(height), diameter(diameter), fluid_holder(type, calc_cylinder(height, diameter), mass)
{
this->level = level;
}

2
src/coolant/condenser.hpp
View File

@ -13,7 +13,7 @@ class condenser : public fluid_holder
public:
condenser(fluid_t type, double height, double diameter, double level);
condenser(fluid_t type, double height, double diameter, double mass, double level);
void update(double dt);
};

37
src/coolant/fluid_holder.cpp
View File

@ -9,7 +9,7 @@
using namespace sim::coolant;
fluid_holder::fluid_holder(fluid_t fluid, double volume) : fluid(fluid), volume(volume)
fluid_holder::fluid_holder(fluid_t fluid, double volume, double extra_mass) : fluid(fluid), volume(volume), extra_mass(extra_mass)
{
}
@ -18,7 +18,7 @@ double fluid_holder::add_heat(double m1, double t1)
{
double t2 = get_heat();
double t = t1 - t2;
double m2 = get_mass() * fluid.jPgk;
double m2 = get_thermal_mass();
double m = m1 + m2;
if(m1 == 0 || m2 == 0)
@ -31,16 +31,18 @@ double fluid_holder::add_heat(double m1, double t1)
double fluid_holder::add_fluid(double v2, double t2)
{
double v1 = level;
double t1 = get_heat();
double t = t1 - t2;
if(level + v2 > volume)
{
v2 = volume - level;
}
heat = t1 - t * v2 / (v1 + v2);
int m1 = get_thermal_mass();
int m2 = fluid.l_to_g(v2);
double t1 = get_heat();
double t = t1 - t2;
heat = t1 - t * m2 / (m1 + m2);
level += v2;
return v2;
@ -64,16 +66,16 @@ double fluid_holder::extract_fluid(double amount)
void fluid_holder::add_steam(double m2, double t2)
{
double m1 = steam;
double m1 = get_thermal_mass();
double t1 = heat;
double m = m1 + m2;
if(m > 0)
{
heat = (t1 * m1 + t2 * m2) / m;
heat = t1 - (t1 - t2) * m2 / (m1 + m2);
}
steam = m;
steam += m2;
}
double fluid_holder::get_pressure() const
@ -81,13 +83,24 @@ double fluid_holder::get_pressure() const
double T = conversions::temperature::c_to_k(heat);
double V = get_steam_volume() * 0.001;
double n = fluid.g_to_mol(steam);
if(V == 0)
{
return NAN;
}
return (n * T * constants::R) / V;
}
double fluid_holder::get_steam_density() const
{
double v = get_steam_volume();
return v > 0 ? steam / v : NAN;
}
void fluid_holder::update(double secs)
{
double mass = get_mass();
double mass = get_thermal_mass();
if(mass > 0)
{
@ -116,7 +129,7 @@ void fluid_holder::update(double secs)
steam = s;
level = fluid.g_to_l(l);
heat -= diff * fluid.jPg / mass / fluid.jPgk;
heat -= diff * fluid.jPg / mass;
}
}

35
src/coolant/fluid_holder.hpp
View File

@ -10,31 +10,34 @@ class fluid_holder
{
protected:
const fluid_t fluid;
const double volume; // litres
double level = 0; // litres
double steam = 0; // grams
double heat = 0; // celsius
fluid_holder(fluid_t fluid, double volume);
public:
double add_heat(double m, double t);
double add_fluid(double amount, double heat);
double extract_fluid(double amount);
fluid_holder(fluid_t fluid, double volume, double extra_mass);
const fluid_t fluid;
const double volume; // litres
const double extra_mass; // grams
virtual double add_heat(double m, double t);
virtual double extract_fluid(double amount);
virtual double add_fluid(double amount, double heat);
virtual void add_steam(double amount, double t);
constexpr double get_volume() const { return volume; } // litres
constexpr double get_level() const { return level; } // litres
constexpr double get_heat() const { return heat; } // celsius
constexpr double get_steam() const { return steam; } // grams
constexpr double get_steam_volume() const { return volume - level; } // litres
constexpr double get_steam_density() const { return steam / get_steam_volume(); } // g/L
constexpr double get_mass() const { return fluid.l_to_g(level) + steam; } // grams
virtual double get_volume() const { return volume; } // litres
virtual double get_level() const { return level; } // litres
virtual double get_heat() const { return heat; } // celsius
virtual double get_steam() const { return steam; } // grams
virtual double get_steam_volume() const { return get_volume() - get_level(); } // litres
virtual double get_mass() const { return fluid.l_to_g(get_level()) + get_steam(); } // grams
virtual double get_thermal_mass() const { return get_mass() + extra_mass; } // grams
virtual double get_pressure() const; // pascals
virtual double get_steam_density() const; // g/L
double get_pressure() const; // pascals
void update(double dt);
};

8
src/coolant/fluid_t.hpp
View File

@ -11,14 +11,12 @@ struct fluid_t
const double gPl; // g/L
const double gPmol; // g/mol
const double jPg; // J/g latent heat of vaporisation
const double jPgk; // J/g/K heat capacity
const double bubble_speed; // m/s
const coolant::vapor_pressure vapor_pressure;
constexpr fluid_t(double gPl, double gPmol, double jPg, double jPgk, double bubble_speed, coolant::vapor_pressure vapor_pressure) :
gPl(gPl), gPmol(gPmol),
jPg(jPg), jPgk(jPgk),
constexpr fluid_t(double gPl, double gPmol, double jPg, double bubble_speed, coolant::vapor_pressure vapor_pressure) :
gPl(gPl), gPmol(gPmol), jPg(jPg),
vapor_pressure(vapor_pressure),
bubble_speed(bubble_speed)
{
@ -33,7 +31,7 @@ struct fluid_t
constexpr double l_to_mol(double l) const { return g_to_mol(l_to_g(l)); }
};
constexpr const fluid_t WATER = fluid_t(1000, 18, 2257, 4.1816, 0.3, {8.07131, 1730.63, 233.426});
constexpr const fluid_t WATER = fluid_t(1000, 18, 2257, 4.1816, {8.07131, 1730.63, 233.426});
}

35
src/coolant/pump.cpp
View File

@ -9,12 +9,12 @@ using namespace sim::coolant;
pump::pump(fluid_holder* src, fluid_holder* dst, double mass, double radius, double l_per_rev, double friction) :
src(src), dst(dst), mass(mass), radius(radius), l_per_rev(l_per_rev), friction(friction)
{
power = 1e3;
power = 1e8;
}
double pump::get_flow() const
{
return l_per_rev * get_rpm() / 60;
return src->fluid.l_to_g(l_per_rev * get_rpm() / 60);
}
double pump::get_rpm() const
@ -22,6 +22,26 @@ double pump::get_rpm() const
return velocity / (M_PI * mass * 0.001 * radius * radius);
}
const char* pump::get_state_string()
{
if(!powered)
{
return "Off";
}
if(idling && std::abs(get_flow()) < 1e-3)
{
return "Idle";
}
if(idling)
{
return "Coasting";
}
return "Revving";
}
static double calc_work(double j, double mass)
{
double m = 1;
@ -36,7 +56,10 @@ static double calc_work(double j, double mass)
void pump::update(double dt)
{
velocity += calc_work(dt * power, mass);
if(powered && !idling)
{
velocity += calc_work(dt * power, mass);
}
double src_heat = src->get_heat();
double p_diff_1 = dst->get_pressure() - src->get_pressure();
@ -54,14 +77,14 @@ void pump::update(double dt)
velocity -= calc_work(work, mass);
if(dst->get_level() > 400 || src->get_level() < 50)
if(dst->get_level() > 400 || src->get_level() < 10)
{
power = 0;
idling = true;
}
else
{
power = 1e3;
idling = false;
}
std::cout << "RPM: " << get_rpm() << "\t";

7
src/coolant/pump.hpp
View File

@ -17,16 +17,19 @@ class pump
const double friction; // J/rev
double velocity = 0; // m/s
double power = 0; // W
public:
double power = 0; // watts
bool powered = true;
bool idling = false;
pump(fluid_holder* src, fluid_holder* dst, double mass, double radius, double l_per_rev, double friction);
double get_flow() const; // L/s
double get_flow() const; // g/s
double get_rpm() const; // rev/min
const char* get_state_string();
void update(double dt);
};

65
src/coolant/valve.cpp
View File

@ -1,14 +1,16 @@
#include "valve.hpp"
#include "../conversions/temperature.hpp"
#include "../constants.hpp"
#include <cmath>
#include <iostream>
using namespace sim::coolant;
valve::valve(fluid_holder* src, fluid_holder* dst, double max) : src(src), dst(dst), max(max)
valve::valve(fluid_holder* src, fluid_holder* dst, double state, double max) : src(src), dst(dst), max(max)
{
this->state = state;
}
void valve::add_open_speed(double v)
@ -28,26 +30,53 @@ void valve::update(double dt)
if(state > 1) state = 1;
if(state < 0) state = 0;
double r = max * state; // L
double m = r * dt;//1 - std::pow(0.5, dt * -std::log2(1 - r));
double pressure1 = src->get_pressure();
double pressure2 = dst->get_pressure();
double density1 = src->get_steam_density(); // g/L
double density2 = dst->get_steam_density(); // g/L
double diff = (pressure1 - pressure2) * m;
double temp, mass;
if(diff > 0)
if(src->get_steam_volume() == 0 || dst->get_steam_volume() == 0)
{
temp = src->get_heat();
mass = diff * src->get_steam_density();
flow = 0;
return;
}
else
double pressure1 = src->get_pressure(); // Pa
double pressure2 = dst->get_pressure();
int overshoots = 0;
double m = max * state * dt;
double temp, mass;
for(;;)
{
temp = dst->get_heat();
mass = diff * dst->get_steam_density();
double diff = (pressure1 - pressure2) * m; // L
if(diff > 0)
{
temp = src->get_heat();
mass = std::min(diff * src->get_steam_density(), src->get_steam());
}
else
{
temp = dst->get_heat();
mass = std::min(diff * dst->get_steam_density(), dst->get_steam());
}
fluid_holder fh_src(*src);
fluid_holder fh_dst(*dst);
fh_src.add_steam(-mass, temp);
fh_dst.add_steam(mass, temp);
// if((pressure1 > fh_dst.get_pressure()) == (pressure2 < fh_src.get_pressure()))
{
break;
}
overshoots += 1;
m *= 0.5;
}
if(overshoots > 0)
{
std::cout << "Warning: overshot " << overshoots << " times\n";
}
src->add_steam(-mass, temp);

2
src/coolant/valve.hpp
View File

@ -19,7 +19,7 @@ class valve
public:
valve(fluid_holder* src, fluid_holder* dst, double max);
valve(fluid_holder* src, fluid_holder* dst, double state, double max);
void update(double secs);
void add_open_speed(double v);

27
src/electric/turbine.cpp
View File

@ -14,32 +14,21 @@ constexpr static double calc_cylinder(double h, double d)
return M_PI * r * r * h;
}
turbine::turbine(coolant::fluid_t type, double height, double diameter, double level) :
height(height), diameter(diameter), sim::coolant::fluid_holder(type, calc_cylinder(height, diameter))
turbine::turbine(coolant::fluid_t type, coolant::condenser* condenser, double length, double diameter, double mass) :
length(length), diameter(diameter), condenser(condenser),
sim::coolant::fluid_holder(type, calc_cylinder(length, diameter), mass)
{
this->level = level;
}
void turbine::update(double secs)
{
sim::system& sys = sim::system::active;
auto& o = *sys.turbine.get();
}
((sim::coolant::fluid_holder*)this)->update(secs);
/* if(level + o.level > o.volume)
{
level += o.level - o.volume;
o.level = o.volume;
}
else
{
o.level += level;
level = 0;
}
balance_steam(o);*/
void turbine::add_steam(double amount, double t)
{
condenser->add_steam(amount, t);
}

22
src/electric/turbine.hpp
View File

@ -2,20 +2,38 @@
#pragma once
#include "../coolant/fluid_holder.hpp"
#include "../coolant/condenser.hpp"
namespace sim::electric
{
class turbine : public sim::coolant::fluid_holder
{
const double height;
coolant::condenser* const condenser;
const double length;
const double diameter;
public:
turbine(coolant::fluid_t type, double height, double diameter, double level);
turbine(coolant::fluid_t type, coolant::condenser* condenser, double length, double diameter, double mass);
void update(double dt);
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); }
virtual double add_fluid(double amount, double heat) { return condenser->add_fluid(amount, heat); }
virtual void add_steam(double amount, double t);
virtual double get_volume() const { return condenser->get_volume(); }
virtual double get_level() const { return condenser->get_level(); }
virtual double get_heat() const { return condenser->get_heat(); } // celsius
virtual double get_steam() const { return condenser->get_steam(); } // grams
virtual double get_steam_volume() const { return condenser->get_steam_volume(); } // litres
virtual double get_mass() const { return condenser->get_mass(); } // grams
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_steam_density() const { return condenser->get_steam_density(); } // g/L
};
};

20
src/graphics/monitor/primary_loop.cpp
View File

@ -73,11 +73,8 @@ void primary_loop::init()
ss << "Opened\nFlow\n\n";
ss << "Turbine Inlet Valve\n\n";
ss << "Opened\nFlow\n\n";
ss << "Turbine\n\n";
ss << "Heat\n";
ss << "Steam\n";
ss << "Pressure\n";
ss << "Level\n\n";
ss << "Primary Pump\n\n";
ss << "Power\nSpeed\nFlow\n\n";
ss << "Condenser\n\n";
ss << "Heat\n";
ss << "Steam\n";
@ -100,20 +97,19 @@ void primary_loop::update()
ss << "\n\n";
ss << show( sys.turbine_bypass_valve->get_state() * 100 ) << " %\n";
ss << show( sys.turbine_bypass_valve->get_flow() ) << " g/s\n";
ss << show( sys.turbine_bypass_valve->get_flow() / 1000 ) << " kg/s\n";
ss << "\n\n\n";
ss << show( sys.turbine_inlet_valve->get_state() * 100 ) << " %\n";
ss << show( sys.turbine_inlet_valve->get_flow() ) << " g/s\n";
ss << show( sys.turbine_inlet_valve->get_flow() / 1000 ) << " kg/s\n";
ss << "\n\n\n";
ss << show( sys.turbine->get_heat() ) << " C\n";
ss << show( sys.turbine->get_steam() ) << " g\n";
ss << show( sys.turbine->get_pressure() / 1000 ) << " kPa\n";
ss << sys.turbine->get_level() << " / " << show( sys.turbine->get_volume() ) << " L\n";
ss << sys.primary_pump->get_state_string() << "\n";
ss << show( sys.primary_pump->get_rpm() ) << " r/min\n";
ss << show( sys.primary_pump->get_flow() / 1000 ) << " kg/s\n";
ss << "\n\n\n";
ss << show( sys.condenser->get_heat() ) << " C\n";
ss << show( sys.condenser->get_steam() ) << " g\n";
ss << show( sys.condenser->get_pressure() / 1000 ) << " kPa\n";
ss << sys.condenser->get_level() << " / " << show( sys.condenser->get_volume() ) << " L\n";
ss << show( sys.condenser->get_level() ) << " / " << show( sys.condenser->get_volume() ) << " L\n";
rmesh.load_text(ss.str().c_str(), 0.04);
mesh2.bind();

4
src/reactor/coolant/vessel.cpp
View File

@ -16,8 +16,8 @@ constexpr static double calc_cylinder(double h, double d)
return M_PI * r * r * h;
}
vessel::vessel(sim::coolant::fluid_t fluid, double height, double diameter, double level) :
sim::coolant::fluid_holder(fluid, calc_cylinder(height, diameter)),
vessel::vessel(sim::coolant::fluid_t fluid, double height, double diameter, double mass, double level) :
sim::coolant::fluid_holder(fluid, calc_cylinder(height, diameter), mass),
height(height), diameter(diameter)
{
this->level = level;

2
src/reactor/coolant/vessel.hpp
View File

@ -16,7 +16,7 @@ public:
const double diameter; // meters
double steam_suspended = 0; // grams
vessel(sim::coolant::fluid_t fluid, double height, double diameter, double level);
vessel(sim::coolant::fluid_t fluid, double height, double diameter, double mass, double level);
double get_steam_suspended() const; // grams
double get_void_ratio() const;

25
src/system.cpp
View File

@ -35,14 +35,14 @@ system::system()
" C C C C "
};
vessel = std::make_unique<reactor::coolant::vessel>(sim::coolant::WATER, 8, 10, 300);
vessel = std::make_unique<reactor::coolant::vessel>(sim::coolant::WATER, 8, 10, 6e6, 300);
reactor = std::make_unique<reactor::reactor>(sim::reactor::builder(19, 19, 1.0 / 4.0, 4, reactor::fuel::fuel_rod(0.5), vessel.get(), layout));
condenser = std::make_unique<coolant::condenser>(sim::coolant::WATER, 8, 6, 0);
turbine = std::make_unique<electric::turbine>(sim::coolant::WATER, 6, 3, 0);
condenser = std::make_unique<coolant::condenser>(sim::coolant::WATER, 8, 6, 3e6, 0);
turbine = std::make_unique<electric::turbine>(sim::coolant::WATER, condenser.get(), 6, 3, 2e6);
turbine_inlet_valve = std::make_unique<coolant::valve>(vessel.get(), turbine.get(), 1e-7);
turbine_bypass_valve = std::make_unique<coolant::valve>(vessel.get(), condenser.get(), 1e-7);
core_pump = std::make_unique<coolant::pump>(condenser.get(), vessel.get(), 1e6, 1, 100, 100);
turbine_inlet_valve = std::make_unique<coolant::valve>(vessel.get(), turbine.get(), 0, 1e-2);
turbine_bypass_valve = std::make_unique<coolant::valve>(vessel.get(), condenser.get(), 0, 1e-2);
primary_pump = std::make_unique<coolant::pump>(condenser.get(), vessel.get(), 1e6, 1, 1, 0);
}
system::system(system&& o)
@ -53,18 +53,21 @@ system::system(system&& o)
turbine = std::move(o.turbine);
turbine_bypass_valve = std::move(o.turbine_bypass_valve);
turbine_inlet_valve = std::move(o.turbine_inlet_valve);
core_pump = std::move(o.core_pump);
primary_pump = std::move(o.primary_pump);
}
void system::update(double dt)
{
dt *= speed;
vessel->update(dt);
reactor->update(dt);
condenser->update(dt);
turbine->update(dt);
vessel->update(dt);
turbine_inlet_valve->update(dt);
turbine_bypass_valve->update(dt);
// core_pump->update(dt);
turbine->update(dt);
condenser->update(dt);
primary_pump->update(dt);
}

2
src/system.hpp
View File

@ -22,7 +22,7 @@ struct system
std::unique_ptr<sim::reactor::coolant::vessel> vessel;
std::unique_ptr<sim::coolant::condenser> condenser;
std::unique_ptr<sim::electric::turbine> turbine;
std::unique_ptr<sim::coolant::pump> core_pump;
std::unique_ptr<sim::coolant::pump> primary_pump;
std::unique_ptr<sim::coolant::valve> turbine_bypass_valve;
std::unique_ptr<sim::coolant::valve> turbine_inlet_valve;
sim::graphics::mesh scene;