using namespace std;
using namespace std::tr1;
+//
+// READ THIS FIRST: physics model used by the steam engine
+//
+// Note: everything here uses SI units unless otherwise stated
+//
+// The "tractive effort" is a measure of the power of a steam engine
+// at a given velocity: P(v). Note: the value usually quoted on the
+// Wikipedia entry for trains is the /starting/ tractive effort
+// (i.e. P(0) = Fmax)
+//
+// Tractive effort is at its maximum value and constant up to some
+// speed TRACTIVE_EFFORT_KNEE after which it decreases as 1/x
+// These values are really engine-dependant, but we're simplifying here.
+// P(v) = { Fmax, if v < TRACTIVE_EFFORT_KNEE
+// { (Fmax * TRACTIVE_EFFORT_KNEE) / v, otherwise
+//
+// Resistance on the train is a combination of friction, drag, and some
+// other sources. This is usually approximated by a quadratic:
+// Q(v) = a + bv + cv^2
+// Where v is velocity. The constants a, b, and c are usually determined
+// experimentally (we'll just have to guess). Where Q(v) intersects the
+// tractive effort curve P(v) determines the train's maximum speed.
+//
+// Run models.gnuplot to see an example of these curves. The functions
+// in that file should match the code here!!
+//
+
// Concrete implementation of a steam engine
class Engine : public IRollingStock,
public IController,
// IController interface
void actOn(Action anAction);
private:
+ double tractiveEffort() const;
+ double resistance() const;
+
IModelPtr myModel;
double mySpeed, myMass, myBoilerPressure, myFireTemp;
double myFuelOnFire;
+ double myStatTractiveEffort;
bool isBrakeOn;
static const double MODEL_SCALE;
+ static const double TRACTIVE_EFFORT_KNEE;
};
const double Engine::MODEL_SCALE(0.4);
+const double Engine::TRACTIVE_EFFORT_KNEE(10.0);
Engine::Engine()
- : mySpeed(0.0), myMass(1000.0), myBoilerPressure(1.0),
- myFireTemp(0.0), myFuelOnFire(0.0), isBrakeOn(true)
+ : mySpeed(0.0), myMass(29.0), myBoilerPressure(1.0),
+ myFireTemp(0.0), myFuelOnFire(0.0),
+ myStatTractiveEffort(34.7),
+ isBrakeOn(true)
{
myModel = loadModel("train.obj", MODEL_SCALE);
}
myModel->render();
}
-// Compute the next state of the engine
-void Engine::update()
+// Calculate the current tractive effort
+double Engine::tractiveEffort() const
{
- const double stopSpeed = 0.001;
-
- // Maximum friction
- const double Fmax = (abs(mySpeed) < stopSpeed ? MU_S : MU_K) * myMass;
-
- // Maximum braking force
- const double Bmax = isBrakeOn ? 2000.0 : 0.0;
+ if (mySpeed < TRACTIVE_EFFORT_KNEE)
+ return myStatTractiveEffort;
+ else
+ return (myStatTractiveEffort * TRACTIVE_EFFORT_KNEE) / mySpeed;
+}
- // `forwards' force
- const double drivingForce = myBoilerPressure;
+// Calculate the magnitude of the resistance on the train
+double Engine::resistance() const
+{
+ const double a = 2.0;
+ const double b = 0.02;
+ const double c = 0.0035;
+ return a + b*mySpeed + c*mySpeed*mySpeed;
+}
- // Net force
- double netForce;
- if (Fmax + Bmax > drivingForce && abs(mySpeed) < stopSpeed)
- netForce = 0.0;
- else
- netForce = drivingForce - Fmax - Bmax;
+// Compute the next state of the engine
+void Engine::update()
+{
+ const double P = tractiveEffort();
+ const double Q = resistance();
- const double accel = netForce / myMass;
+ // A fudge factor to make the acceleration look realistic
+ const double MASS_TWEAK = 10.0;
- // Consume some fuel
- const double burnRate = 0.999;
- myFuelOnFire *= burnRate;
-
- // The fire temparature is simply a function of fuel
- // TODO: find a better function!
- myFireTemp = min(myFuelOnFire * 100.0, 1000.0);
-
- // Boiler pressure is roughly proportional to temparature
- myBoilerPressure = myFireTemp * 10.0;
-
- // Accelerate the train
- mySpeed += accel / 1000.0;
-
- // Apply drag: drag is roughly proportional to the square of
- // velocity at high speeds
- const double dragCoeff = 0.1;
- mySpeed -= dragCoeff * mySpeed * mySpeed;
-
- debug() << "Hold: " << (Fmax + Bmax)
- << ", go: " << drivingForce
- << ", temp=" << myFireTemp
- << ", pressure=" << myBoilerPressure
- << ", accel=" << accel
- << ", drag=" << (dragCoeff * mySpeed * mySpeed);
+ const double a = (P - Q) / (myMass * MASS_TWEAK);
+
+ mySpeed += a;
+ debug() << "P=" << P << ", Q=" << Q
+ << ", a=" << a << ", v=" << mySpeed;
}
// User interface to the engine