_vstabs.add(vstab);
}
-void Airplane::addFuselage(float* front, float* back, float width)
+void Airplane::addFuselage(float* front, float* back, float width,
+ float taper, float mid)
{
Fuselage* f = new Fuselage();
int i;
f->back[i] = back[i];
}
f->width = width;
+ f->taper = taper;
+ f->mid = mid;
_fuselages.add(f);
}
int i;
for(i=0; i<w->numSurfaces(); i++) {
Surface* s = (Surface*)w->getSurface(i);
+
+ float td = s->getTotalDrag();
+ s->setTotalDrag(td);
+
_model.addSurface(s);
+ float mass = w->getSurfaceWeight(i);
+ mass = mass * Math::sqrt(mass);
float pos[3];
s->getPosition(pos);
- _model.getBody()->addMass(w->getSurfaceWeight(i), pos);
- wgt += w->getSurfaceWeight(i);
+ _model.getBody()->addMass(mass, pos);
+ wgt += mass;
}
return wgt;
}
int j;
for(j=0; j<segs; j++) {
float frac = (j+0.5) / segs;
+
+ float scale = 1;
+ if(frac < f->mid)
+ scale = f->taper+(1-f->taper) * (frac / f->mid);
+ else
+ scale = f->taper+(1-f->taper) * (frac - f->mid) / (1 - f->mid);
+
+ // Where are we?
float pos[3];
Math::mul3(frac, fwd, pos);
Math::add3(f->back, pos, pos);
- _model.getBody()->addMass(segWgt, pos);
- wgt += segWgt;
+
+ // _Mass_ weighting goes as surface area^(3/2)
+ float mass = scale*segWgt * Math::sqrt(scale*segWgt);
+ _model.getBody()->addMass(mass, pos);
+ wgt += mass;
// Make a Surface too
Surface* s = new Surface();
float sideDrag = len/wid;
s->setYDrag(sideDrag);
s->setZDrag(sideDrag);
- s->setTotalDrag(segWgt);
+ s->setTotalDrag(scale*segWgt);
// FIXME: fails for fuselages aligned along the Y axis
float o[9];
void setTail(Wing* tail);
void addVStab(Wing* vstab);
- void addFuselage(float* front, float* back, float width);
+ void addFuselage(float* front, float* back, float width,
+ float taper=1, float mid=0.5);
int addTank(float* pos, float cap, float fuelDensity);
void addGear(Gear* g, float transitionTime);
void addThruster(Thruster* t, float mass, float* cg);
float getFuelDensity(int tank); // kg/m^3
void compile(); // generate point masses & such, then solve
+ void stabilizeThrust();
// Solution output values
int getSolutionIterations();
private:
struct Tank { float pos[3]; float cap; float fill;
float density; int handle; };
- struct Fuselage { float front[3], back[3], width; };
+ struct Fuselage { float front[3], back[3], width, taper, mid; };
struct GearRec { Gear* gear; Surface* surf; float wgt; float time; };
struct ThrustRec { Thruster* thruster;
int handle; float cg[3]; float mass; };
float compileWing(Wing* w);
float compileFuselage(Fuselage* f);
void compileGear(GearRec* gr);
- void stabilizeThrust();
void applyDragFactor(float factor);
void applyLiftRatio(float factor);
float clamp(float val, float min, float max);
// P in pascals (N/m^2), rho is kg/m^3, T in kelvin.
const float R = 287.1;
+// Specific heat ratio for air, at "low" temperatures.
+const float GAMMA = 1.4;
+
float Atmosphere::getStdTemperature(float alt)
{
return getRecord(alt, 1);
float Atmosphere::calcMach(float spd, float temp)
{
- return spd / Math::sqrt(1.4 * R * temp);
+ return spd / Math::sqrt(GAMMA * R * temp);
+}
+
+void Atmosphere::calcStaticAir(float p0, float t0, float d0, float v,
+ float* pOut, float* tOut, float* dOut)
+{
+ const static float C0 = ((GAMMA-1)/(2*R*GAMMA));
+ const static float C1 = 1/(GAMMA-1);
+
+ *tOut = t0 + (v*v) * C0;
+ *dOut = d0 * Math::pow(*tOut / t0, C1);
+ *pOut = (*dOut) * R * (*tOut);
}
float Atmosphere::getRecord(float alt, int recNum)
static float calcVEAS(float spd, float pressure, float temp);
static float calcMach(float spd, float temp);
static float calcDensity(float pressure, float temp);
+
+ // Given ambient ("0") pressure/density/temperature values,
+ // calculate the properties of static air (air accelerated to the
+ // aircraft's speed) at a given velocity. Includes
+ // compressibility, but not shock effects.
+ static void calcStaticAir(float p0, float t0, float d0, float v,
+ float* pOut, float* tOut, float* dOut);
private:
static float getRecord(float alt, int idx);
#include "Jet.hpp"
#include "Thruster.hpp"
#include "PropEngine.hpp"
+#include "PistonEngine.hpp"
#include "Gear.hpp"
#include "Wing.hpp"
#include "Math.hpp"
delete v;
}
- for(i=0; i<_outputs.size(); i++) {
- OutRec* o = (OutRec*)_outputs.get(i);
- delete[] o->values;
- delete o;
- }
+ for(i=0; i<_outputs.size(); i++)
+ delete (OutRec*)_outputs.get(i);
}
int ControlMap::newInput()
return _inputs.add(v);
}
+void ControlMap::addMapping(int input, int type, void* object, int options,
+ float src0, float src1, float dst0, float dst1)
+{
+ addMapping(input, type, object, options);
+
+ // The one we just added is last in the list (ugly, awful hack!)
+ Vector* maps = (Vector*)_inputs.get(input);
+ MapRec* m = (MapRec*)maps->get(maps->size() - 1);
+
+ m->src0 = src0;
+ m->src1 = src1;
+ m->dst0 = dst0;
+ m->dst1 = dst1;
+}
+
void ControlMap::addMapping(int input, int type, void* object, int options)
{
// See if the output object already exists
out = new OutRec();
out->type = type;
out->object = object;
- out->n = 0;
- out->values = 0;
_outputs.add(out);
}
-
- // Make space for the new input value
- int idx = out->n++;
- delete[] out->values;
- out->values = new float[out->n];
-
- // Add the new option tag
- out->options.add((void*)options);
-
+
// Make a new input record
MapRec* map = new MapRec();
map->out = out;
- map->idx = idx;
+ map->opt = options;
+ map->idx = out->maps.add(map);
+
+ // The default ranges differ depending on type!
+ map->src1 = map->dst1 = 1;
+ map->src0 = map->dst0 = 0;
+ if(type==FLAP0 || type==FLAP1 || type==STEER)
+ map->src0 = map->dst0 = -1;
- // And add it to the approproate vector.
+ // And add it to the approproate vectors.
Vector* maps = (Vector*)_inputs.get(input);
maps->add(map);
}
void ControlMap::reset()
{
// Set all the values to zero
- int i;
- for(i=0; i<_outputs.size(); i++) {
+ for(int i=0; i<_outputs.size(); i++) {
OutRec* o = (OutRec*)_outputs.get(i);
- int j;
- for(j=0; j<o->n; j++)
- o->values[j] = 0;
+ for(int j=0; j<o->maps.size(); j++)
+ ((MapRec*)o->maps.get(j))->val = 0;
}
}
-void ControlMap::setInput(int input, float value)
+void ControlMap::setInput(int input, float val)
{
Vector* maps = (Vector*)_inputs.get(input);
- int i;
- for(i=0; i<maps->size(); i++) {
- MapRec* map = (MapRec*)maps->get(i);
- map->out->values[map->idx] = value;
+ for(int i=0; i<maps->size(); i++) {
+ MapRec* m = (MapRec*)maps->get(i);
+
+ float val2 = val;
+
+ // Do the scaling operation. Clamp to [src0:src1], rescale to
+ // [0:1] within that range, then map to [dst0:dst1].
+ if(val2 < m->src0) val2 = m->src0;
+ if(val2 > m->src1) val2 = m->src1;
+ val2 = (val2 - m->src0) / (m->src1 - m->src0);
+ val2 = m->dst0 + val2 * (m->dst1 - m->dst0);
+
+ m->val = val2;
}
}
// control axes like ailerons.
float lval = 0, rval = 0;
int i;
- for(i=0; i<o->n; i++) {
- float val = o->values[i];
- int opt = (int)o->options.get(i);
- if(opt & OPT_SQUARE)
+ for(i=0; i<o->maps.size(); i++) {
+ MapRec* m = (MapRec*)o->maps.get(i);
+ float val = m->val;
+
+ if(m->opt & OPT_SQUARE)
val = val * Math::abs(val);
- if(opt & OPT_INVERT)
+ if(m->opt & OPT_INVERT)
val = -val;
lval += val;
- if(opt & OPT_SPLIT)
+ if(m->opt & OPT_SPLIT)
rval -= val;
else
rval += val;
case MIXTURE: ((Thruster*)obj)->setMixture(lval); break;
case ADVANCE: ((PropEngine*)obj)->setAdvance(lval); break;
case REHEAT: ((Jet*)obj)->setReheat(lval); break;
+ case VECTOR: ((Jet*)obj)->setRotation(lval); break;
case BRAKE: ((Gear*)obj)->setBrake(lval); break;
case STEER: ((Gear*)obj)->setRotation(lval); break;
case EXTEND: ((Gear*)obj)->setExtension(lval); break;
case FLAP0: ((Wing*)obj)->setFlap0(lval, rval); break;
case FLAP1: ((Wing*)obj)->setFlap1(lval, rval); break;
case SPOILER: ((Wing*)obj)->setSpoiler(lval, rval); break;
+ case BOOST:
+ ((Thruster*)obj)->getPistonEngine()->setBoost(lval);
+ break;
}
}
}
enum OutputType { THROTTLE, MIXTURE, ADVANCE, REHEAT, PROP,
BRAKE, STEER, EXTEND,
- INCIDENCE, FLAP0, FLAP1, SLAT, SPOILER };
+ INCIDENCE, FLAP0, FLAP1, SLAT, SPOILER, VECTOR,
+ BOOST };
enum { OPT_SPLIT = 0x01,
OPT_INVERT = 0x02,
// of object!
void addMapping(int input, int output, void* object, int options=0);
+ // An additional form to specify a mapping range. Input values
+ // outside of [src0:src1] are clamped, and are then mapped to
+ // [dst0:dst1] before being set on the object.
+ void addMapping(int input, int output, void* object, int options,
+ float src0, float src1, float dst0, float dst1);
+
// Resets our accumulated input values. Call before any
// setInput() invokations.
void reset();
void applyControls();
private:
- struct OutRec { int type; void* object; int n;
- float* values; Vector options; };
- struct MapRec { OutRec* out; int idx; };
+ struct OutRec { int type; void* object; Vector maps; };
+ struct MapRec { OutRec* out; int idx; int opt; float val;
+ float src0; float src1; float dst0; float dst1; };
// A list of (sub)Vectors containing a bunch of MapRec objects for
// each input handle.
#include <Main/fg_props.hxx>
#include "Jet.hpp"
+#include "SimpleJet.hpp"
#include "Gear.hpp"
#include "Atmosphere.hpp"
#include "PropEngine.hpp"
static const float CM2GALS = 264.172037284;
static const float HP2W = 745.700;
static const float INHG2PA = 3386.389;
+static const float K2DEGF = 1.8;
+static const float CIN2CM = 1.6387064e-5;
// Stubs, so that this can be compiled without the FlightGear
// binary. What's the best way to handle this?
delete[] a->name;
delete a;
}
- for(i=0; i<_pistons.size(); i++) {
- EngRec* er = (EngRec*)_pistons.get(i);
+ for(i=0; i<_thrusters.size(); i++) {
+ EngRec* er = (EngRec*)_thrusters.get(i);
delete[] er->prefix;
- delete (PropEngine*)er->eng;
- delete er;
- }
- for(i=0; i<_jets.size(); i++) {
- EngRec* er = (EngRec*)_pistons.get(i);
- delete[] er->prefix;
- delete (Jet*)er->eng;
+ delete er->eng;
delete er;
}
for(i=0; i<_weights.size(); i++) {
void FGFDM::init()
{
- // We don't want to use these ties (we set the values ourselves)
- fgUntie("/consumables/fuel/tank[0]/level-gal_us");
- fgUntie("/consumables/fuel/tank[1]/level-gal_us");
-
// Allows the user to start with something other than full fuel
- _airplane.setFuelFraction(fgGetFloat("/yasim/fuel-fraction", 1));
+ _airplane.setFuelFraction(fgGetFloat("/sim/fuel-fraction", 1));
// This has a nasty habit of being false at startup. That's not
// good.
_airplane.addVStab(parseWing(a, name));
} else if(eq(name, "propeller")) {
parsePropeller(a);
+ } else if(eq(name, "thruster")) {
+ SimpleJet* j = new SimpleJet();
+ _currObj = j;
+ v[0] = attrf(a, "x"); v[1] = attrf(a, "y"); v[2] = attrf(a, "z");
+ j->setPosition(v);
+ _airplane.addThruster(j, 0, v);
+ v[0] = attrf(a, "vx"); v[1] = attrf(a, "vy"); v[2] = attrf(a, "vz");
+ j->setDirection(v);
+ j->setThrust(attrf(a, "thrust") * LBS2N);
} else if(eq(name, "jet")) {
Jet* j = new Jet();
_currObj = j;
v[1] = attrf(a, "y");
v[2] = attrf(a, "z");
float mass = attrf(a, "mass") * LBS2KG;
- j->setDryThrust(attrf(a, "thrust") * LBS2N);
- j->setReheatThrust(attrf(a, "afterburner", 0) * LBS2N);
+ j->setMaxThrust(attrf(a, "thrust") * LBS2N,
+ attrf(a, "afterburner", 0) * LBS2N);
+ j->setVectorAngle(attrf(a, "rotate", 0) * DEG2RAD);
+
+ float n1min = attrf(a, "n1-idle", 55);
+ float n1max = attrf(a, "n1-max", 102);
+ float n2min = attrf(a, "n2-idle", 73);
+ float n2max = attrf(a, "n2-max", 103);
+ j->setRPMs(n1min, n1max, n2min, n2max);
+
+ if(a->hasAttribute("tsfc")) j->setTSFC(attrf(a, "tsfc"));
+ if(a->hasAttribute("egt")) j->setEGT(attrf(a, "egt"));
+ if(a->hasAttribute("epr")) j->setEPR(attrf(a, "epr"));
+ if(a->hasAttribute("exhaust-speed"))
+ j->setVMax(attrf(a, "exhaust-speed") * KTS2MPS);
+
j->setPosition(v);
_airplane.addThruster(j, mass, v);
sprintf(buf, "/engines/engine[%d]", _nextEngine++);
EngRec* er = new EngRec();
er->eng = j;
er->prefix = dup(buf);
- _jets.add(er);
+ _thrusters.add(er);
} else if(eq(name, "gear")) {
Gear* g = new Gear();
_currObj = g;
g->setCompression(v);
g->setStaticFriction(attrf(a, "sfric", 0.8));
g->setDynamicFriction(attrf(a, "dfric", 0.7));
+ if(a->hasAttribute("castering"))
+ g->setCastering(true);
float transitionTime = attrf(a, "retract-time", 0);
_airplane.addGear(g, transitionTime);
} else if(eq(name, "fuselage")) {
b[0] = attrf(a, "bx");
b[1] = attrf(a, "by");
b[2] = attrf(a, "bz");
- _airplane.addFuselage(v, b, attrf(a, "width"));
+ float taper = attrf(a, "taper", 1);
+ float mid = attrf(a, "midpoint", 0.5);
+ _airplane.addFuselage(v, b, attrf(a, "width"), taper, mid);
} else if(eq(name, "tank")) {
v[0] = attrf(a, "x");
v[1] = attrf(a, "y");
v[2] = attrf(a, "z");
float density = 6.0; // gasoline, in lbs/gal
if(a->hasAttribute("jet")) density = 6.72;
- density *= LBS2KG/CM2GALS;
+ density *= LBS2KG*CM2GALS;
_airplane.addTank(v, attrf(a, "capacity") * LBS2KG, density);
} else if(eq(name, "ballast")) {
v[0] = attrf(a, "x");
opt |= a->hasAttribute("split") ? ControlMap::OPT_SPLIT : 0;
opt |= a->hasAttribute("invert") ? ControlMap::OPT_INVERT : 0;
opt |= a->hasAttribute("square") ? ControlMap::OPT_SQUARE : 0;
- _airplane.getControlMap()->addMapping(parseAxis(axis),
- parseOutput(output),
- _currObj,
- opt);
+
+ ControlMap* cm = _airplane.getControlMap();
+ if(a->hasAttribute("src0")) {
+ cm->addMapping(parseAxis(axis), parseOutput(output),
+ _currObj, opt,
+ attrf(a, "src0"), attrf(a, "src1"),
+ attrf(a, "dst0"), attrf(a, "dst1"));
+ } else {
+ cm->addMapping(parseAxis(axis), parseOutput(output),
+ _currObj, opt);
+ }
} else {
// assert: must be under a "cruise" or "approach" tag
float value = attrf(a, "value", 0);
{
char buf[256];
int i;
+ float fuelDensity = 718.95; // default to gasoline: ~6 lb/gal
for(i=0; i<_airplane.numTanks(); i++) {
+ fuelDensity = _airplane.getFuelDensity(i);
sprintf(buf, "/consumables/fuel/tank[%d]/level-gal_us", i);
- fgSetFloat(buf,
- CM2GALS*_airplane.getFuel(i)/_airplane.getFuelDensity(i));
+ fgSetFloat(buf, CM2GALS*_airplane.getFuel(i)/fuelDensity);
}
- for(i=0; i<_pistons.size(); i++) {
- EngRec* er = (EngRec*)_pistons.get(i);
- PropEngine* p = (PropEngine*)er->eng;
-
- sprintf(buf, "%s/rpm", er->prefix);
- fgSetFloat(buf, p->getOmega() / RPM2RAD);
+ for(i=0; i<_thrusters.size(); i++) {
+ EngRec* er = (EngRec*)_thrusters.get(i);
+ Thruster* t = er->eng;
sprintf(buf, "%s/fuel-flow-gph", er->prefix);
- fgSetFloat(buf, p->getFuelFlow() * (3600*2.2/5)); // FIXME, wrong
- }
+ fgSetFloat(buf, (t->getFuelFlow()/fuelDensity) * 3600 * CM2GALS);
- for(i=0; i<_jets.size(); i++) {
- EngRec* er = (EngRec*)_jets.get(i);
- Jet* j = (Jet*)er->eng;
-
- sprintf(buf, "%s/fuel-flow-gph", er->prefix);
- fgSetFloat(buf, j->getFuelFlow() * (3600*2.2/6)); // FIXME, wrong
+ if(t->getPropEngine()) {
+ PropEngine* p = t->getPropEngine();
+
+ sprintf(buf, "%s/rpm", er->prefix);
+ fgSetFloat(buf, p->getOmega() / RPM2RAD);
+ }
+
+ if(t->getPistonEngine()) {
+ PistonEngine* p = t->getPistonEngine();
+
+ sprintf(buf, "%s/mp-osi", er->prefix);
+ fgSetFloat(buf, p->getMP() * (1/INHG2PA));
+
+ sprintf(buf, "%s/egt-degf", er->prefix);
+ fgSetFloat(buf, p->getEGT() * K2DEGF + 459.4);
+ }
+
+ if(t->getJet()) {
+ Jet* j = t->getJet();
+
+ sprintf(buf, "%s/n1", er->prefix);
+ fgSetFloat(buf, j->getN1());
+
+ sprintf(buf, "%s/n2", er->prefix);
+ fgSetFloat(buf, j->getN2());
+
+ sprintf(buf, "%s/epr", er->prefix);
+ fgSetFloat(buf, j->getEPR());
+
+ sprintf(buf, "%s/egt-degf", er->prefix);
+ fgSetFloat(buf, j->getEGT() * K2DEGF + 459.4);
+ }
}
}
PropEngine* thruster = new PropEngine(prop, eng, moment);
_airplane.addThruster(thruster, mass, cg);
+ if(a->hasAttribute("displacement"))
+ eng->setDisplacement(attrf(a, "displacement") * CIN2CM);
+
+ if(a->hasAttribute("compression"))
+ eng->setCompression(attrf(a, "compression"));
+
if(a->hasAttribute("turbo-mul")) {
float mul = attrf(a, "turbo-mul");
float mp = attrf(a, "wastegate-mp", 1e6) * INHG2PA;
EngRec* er = new EngRec();
er->eng = thruster;
er->prefix = dup(buf);
- _pistons.add(er);
+ _thrusters.add(er);
_currObj = thruster;
}
if(eq(name, "MIXTURE")) return ControlMap::MIXTURE;
if(eq(name, "ADVANCE")) return ControlMap::ADVANCE;
if(eq(name, "REHEAT")) return ControlMap::REHEAT;
+ if(eq(name, "BOOST")) return ControlMap::BOOST;
+ if(eq(name, "VECTOR")) return ControlMap::VECTOR;
if(eq(name, "PROP")) return ControlMap::PROP;
if(eq(name, "BRAKE")) return ControlMap::BRAKE;
if(eq(name, "STEER")) return ControlMap::STEER;
if(eq(name, "SLAT")) return ControlMap::SLAT;
if(eq(name, "SPOILER")) return ControlMap::SPOILER;
// error here...
- return -1;
+ return *(int*)0;
}
void FGFDM::parseWeight(XMLAttributes* a)
~FGFDM();
void init();
void iterate(float dt);
+ void getExternalInput(float dt=1e6);
Airplane* getAirplane();
private:
struct AxisRec { char* name; int handle; };
- struct EngRec { char* prefix; void* eng; };
+ struct EngRec { char* prefix; Thruster* eng; };
struct WeightRec { char* prop; float size; int handle; };
- void getExternalInput(float dt);
void setOutputProperties();
Wing* parseWing(XMLAttributes* a, const char* name);
Vector _weights;
// Engine types. Contains an EngRec structure.
- Vector _jets;
- Vector _pistons;
+ Vector _thrusters;
// Parsing temporaries
void* _currObj;
_brake = 0;
_rot = 0;
_extension = 1;
+ _castering = false;
}
void Gear::setPosition(float* position)
_extension = Math::clamp(extension, 0, 1);
}
+void Gear::setCastering(bool c)
+{
+ _castering = c;
+}
+
void Gear::getPosition(float* out)
{
int i;
return _frac;
}
+bool Gear::getCastering()
+{
+ return _castering;
+}
+
void Gear::calcForce(RigidBody* body, float* v, float* rot, float* ground)
{
// Init the return values
_wow = (fmag - damp) * -Math::dot3(cmpr, ground);
Math::mul3(-_wow, ground, _force);
+ // Castering gear feel no force in the ground plane
+ if(_castering)
+ return;
+
// Wheels are funky. Split the velocity along the ground plane
// into rolling and skidding components. Assuming small angles,
// we generate "forward" and "left" unit vectors (the compression
float Gear::calcFriction(float wgt, float v)
{
- // How slow is stopped? 50 cm/second?
- const float STOP = 0.5;
+ // How slow is stopped? 10 cm/second?
+ const float STOP = 0.1;
const float iSTOP = 1/STOP;
v = Math::abs(v);
if(v < STOP) return v*iSTOP * wgt * _sfric;
void setBrake(float brake);
void setRotation(float rotation);
void setExtension(float extension);
-
+ void setCastering(bool castering);
+
void getPosition(float* out);
void getCompression(float* out);
float getSpring();
float getBrake();
float getRotation();
float getExtension();
+ bool getCastering();
// Takes a velocity of the aircraft relative to ground, a rotation
// vector, and a ground plane (all specified in local coordinates)
private:
float calcFriction(float wgt, float v);
+ bool _castering;
float _pos[3];
float _cmpr[3];
float _spring;
Jet::Jet()
{
- _rho0 = Atmosphere::getStdDensity(0);
- _thrust = 0;
- _abThrust = 0;
+ _maxThrust = 0;
+ _abFactor = 1;
_reheat = 0;
+ _rotControl = 0;
+ _maxRot = 0;
+
+ // Initialize parameters for an early-ish subsonic turbojet. More
+ // recent turbofans will typically have a lower vMax, epr0, and
+ // tsfc.
+ _vMax = 800;
+ _epr0 = 3.0;
+ _tsfc = 0.8;
+ _egt0 = 1050;
+ _n1Min = 55;
+ _n1Max = 102;
+ _n2Min = 73;
+ _n2Max = 103;
+ setSpooling(4); // 4 second spool time? s'bout right.
+
+ // And initialize to an engine that is idling
+ _n1 = _n1Min;
+ _n2 = _n2Min;
+
+ // And sanify the remaining junk, just in case.
+ _epr = 1;
+ _fuelFlow = 0;
+ _egt = 273;
+ _tempCorrect = 1;
+ _pressureCorrect = 1;
}
void Jet::stabilize()
{
- return; // no-op for now
+ // Just run it for an hour, there's no need to iterate given the
+ // algorithms used.
+ integrate(3600);
+}
+
+void Jet::setMaxThrust(float thrust, float afterburner)
+{
+ _maxThrust = thrust;
+ if(afterburner == 0) _abFactor = 1;
+ else _abFactor = afterburner/thrust;
+}
+
+void Jet::setVMax(float spd)
+{
+ _vMax = spd;
+}
+
+void Jet::setTSFC(float tsfc)
+{
+ _tsfc = tsfc;
}
-void Jet::setDryThrust(float thrust)
+void Jet::setRPMs(float idleN1, float maxN1, float idleN2, float maxN2)
{
- _thrust = thrust;
+ _n1Min = idleN1;
+ _n1Max = maxN1;
+ _n2Min = idleN2;
+ _n2Max = maxN2;
}
-void Jet::setReheatThrust(float thrust)
+void Jet::setEGT(float takeoffEGT)
{
- _abThrust = thrust;
+ _egt0 = takeoffEGT;
+}
+
+void Jet::setEPR(float takeoffEPR)
+{
+ _epr0 = takeoffEPR;
+}
+
+void Jet::setSpooling(float time)
+{
+ // 2.3 = -ln(0.1), i.e. x=2.3 is the 90% point we're defining
+ // The extra fudge factor is there because the N1 speed (which
+ // determines thrust) lags the N2 speed.
+ _decay = 1.5 * 2.3 / time;
+}
+
+void Jet::setVectorAngle(float angle)
+{
+ _maxRot = angle;
}
void Jet::setReheat(float reheat)
_reheat = Math::clamp(reheat, 0, 1);
}
-void Jet::getThrust(float* out)
+void Jet::setRotation(float rot)
{
- float t = _thrust * _throttle;
- t += (_abThrust - _thrust) * _reheat;
- t *= _rho / _rho0;
- Math::mul3(t, _dir, out);
+ if(rot < 0) rot = 0;
+ if(rot > 1) rot = 1;
+ _rotControl = rot;
}
-void Jet::getTorque(float* out)
+
+float Jet::getN1()
{
- out[0] = out[1] = out[2] = 0;
- return;
+ return _n1 * _tempCorrect;
}
-void Jet::getGyro(float* out)
+float Jet::getN2()
{
- out[0] = out[1] = out[2] = 0;
- return;
+ return _n2 * _tempCorrect;
+}
+
+float Jet::getEPR()
+{
+ return _epr;
+}
+
+float Jet::getEGT()
+{
+ // Exactly zero means "off" -- return the ambient temperature
+ if(_egt == 0) return _temp;
+
+ return _egt * _tempCorrect * _tempCorrect;
}
float Jet::getFuelFlow()
{
- return 0;
+ return _fuelFlow * _pressureCorrect;
}
void Jet::integrate(float dt)
{
+ // Sea-level values
+ const static float P0 = Atmosphere::getStdPressure(0);
+ const static float T0 = Atmosphere::getStdTemperature(0);
+ const static float D0 = Atmosphere::getStdDensity(0);
+
+ float speed = -Math::dot3(_wind, _dir);
+
+ float statT, statP, statD;
+ Atmosphere::calcStaticAir(_pressure, _temp, _rho, speed,
+ &statP, &statT, &statD);
+ _pressureCorrect = statP/P0;
+ _tempCorrect = Math::sqrt(statT/T0);
+
+ // Linearly taper maxThrust to zero at vMax
+ float vCorr = 1 - (speed/_vMax);
+
+ float maxThrust = _maxThrust * vCorr * (statD/D0);
+ _thrust = maxThrust * _throttle;
+
+ // Now get a "beta" (i.e. EPR - 1) value. The output values are
+ // expressed as functions of beta.
+ float ibeta0 = 1/(_epr0 - 1);
+ float betaTarget = (_epr0 - 1) * (_thrust/_maxThrust) * (P0/_pressure)
+ * (_temp/statT);
+ float n2Target = _n2Min + (betaTarget*ibeta0) * (_n2Max - _n2Min);
+
+ // Note that this "first" beta value is used to compute a target
+ // for N2 only Integrate the N2 speed and back-calculate a beta1
+ // target. The N1 speed will seek to this.
+ _n2 = (_n2 + dt*_decay * n2Target) / (1 + dt*_decay);
+
+ float betaN2 = (_epr0-1) * (_n2 - _n2Min) / (_n2Max - _n2Min);
+ float n1Target = _n1Min + betaN2*ibeta0 * (_n1Max - _n1Min);
+ _n1 = (_n1 + dt*_decay * n1Target) / (1 + dt*_decay);
+
+ // The actual thrust produced is keyed to the N1 speed. Add the
+ // afterburners in at the end.
+ float betaN1 = (_epr0-1) * (_n1 - _n1Min) / (_n1Max - _n1Min);
+ _thrust *= betaN1/(betaTarget+.00001); // blowup protection
+ _thrust *= 1 + _reheat*(_abFactor-1);
+
+ // Finally, calculate the output variables. Use a 80/20 mix of
+ // the N2/N1 speeds as the key.
+ float beta = 0.8*betaN2 + 0.2*betaN1;
+ _epr = beta + 1;
+ float ff0 = _maxThrust*_tsfc*(1/(3600*9.8)); // takeoff fuel flow, kg/s
+ _fuelFlow = ff0 * beta*ibeta0;
+ _fuelFlow *= 1 + (3.5 * _reheat * _abFactor); // Afterburners take
+ // 3.5 times as much
+ // fuel per thrust unit
+ _egt = T0 + beta*ibeta0 * (_egt0 - T0);
+}
+
+void Jet::getThrust(float* out)
+{
+ Math::mul3(_thrust, _dir, out);
+
+ // Rotate about the Y axis for thrust vectoring
+ float angle = _rotControl * _maxRot;
+ float s = Math::sin(angle);
+ float c = Math::cos(angle);
+ float o0 = out[0];
+ out[0] = c * o0 + s * out[2];
+ out[2] = -s * o0 + c * out[2];
+}
+
+void Jet::getTorque(float* out)
+{
+ out[0] = out[1] = out[2] = 0;
+ return;
+}
+
+void Jet::getGyro(float* out)
+{
+ out[0] = out[1] = out[2] = 0;
return;
}
namespace yasim {
-// Incredibly dumb placeholder for a Jet implementation. But it does
-// what's important, which is provide thrust.
class Jet : public Thruster {
public:
Jet();
virtual Jet* getJet() { return this; }
+
+ void setMaxThrust(float thrust, float afterburner=0);
+ void setVMax(float spd);
+ void setTSFC(float tsfc);
+ void setRPMs(float idleN1, float maxN1, float idleN2, float maxN2);
+ void setEGT(float takeoffEGT);
+ void setEPR(float takeoffEPR);
+ void setVectorAngle(float angle);
- void setDryThrust(float thrust);
- void setReheatThrust(float thrust);
+ // The time it takes the engine to reach 90% thrust from idle
+ void setSpooling(float time);
+
+ // Sets the reheat control
void setReheat(float reheat);
+ // Sets the thrust vector control (0-1)
+ void setRotation(float rot);
+
+ float getN1();
+ float getN2();
+ float getEPR();
+ float getEGT();
+
+ // From Thruster:
virtual void getThrust(float* out);
virtual void getTorque(float* out);
virtual void getGyro(float* out);
virtual void stabilize();
private:
- float _thrust;
- float _abThrust;
- float _rho0;
float _reheat;
+
+ float _maxThrust; // Max dry thrust at sea level
+ float _abFactor; // Afterburner thrust multiplier
+
+ float _maxRot;
+ float _rotControl;
+
+ float _decay; // time constant for the exponential integration
+ float _vMax; // speed at which thrust is zero
+ float _epr0; // EPR at takeoff thrust
+ float _tsfc; // TSFC ((lb/hr)/lb) at takeoff thrust and zero airspeed
+ float _egt0; // EGT at takeoff thrust
+ float _n1Min; // N1 at ground idle
+ float _n1Max; // N1 at takeoff thrust
+ float _n2Min; // N2 at ground idle
+ float _n2Max; // N2 at takeoff thrust
+
+ float _thrust; // Current thrust
+ float _epr; // Current EPR
+ float _n1; // Current UNCORRECTED N1 (percent)
+ float _n2; // Current UNCORRECTED N2 (percent)
+ float _fuelFlow; // Current UNCORRECTED fuel flow (kg/s)
+ float _egt; // Current UNCORRECTED EGT (kelvin)
+
+ float _tempCorrect; // Intake temp / std temp (273 K)
+ float _pressureCorrect; // Intake pressure / std pressure
};
}; // namespace yasim
FGFDM.cpp Gear.cpp Glue.cpp Integrator.cpp Jet.cpp \
Math.cpp Model.cpp PistonEngine.cpp Propeller.cpp \
PropEngine.cpp RigidBody.cpp Surface.cpp \
- Thruster.cpp Wing.cpp
+ Thruster.cpp Wing.cpp SimpleJet.cpp
INCLUDES += -I$(top_srcdir) -I$(top_srcdir)/src
#include "PistonEngine.hpp"
namespace yasim {
+const static float HP2W = 745.7;
+const static float CIN2CM = 1.6387064e-5;
+
PistonEngine::PistonEngine(float power, float speed)
{
+ _boost = 1;
+
// Presume a BSFC (in lb/hour per HP) of 0.45. In SI that becomes
- // (2.2 lb/kg, 745.7 W/hp, 3600 sec/hour) 3.69e-07 kg/Ws.
- _f0 = power * 3.69e-07;
+ // (2.2 lb/kg, 745.7 W/hp, 3600 sec/hour) 7.62e-08 kg/Ws.
+ _f0 = power * 7.62e-08;
_power0 = power;
_omega0 = speed;
_turbo = 1;
_maxMP = 1e6; // No waste gate on non-turbo engines.
+
+ // Guess at reasonable values for these guys. Displacements run
+ // at about 2 cubic inches per horsepower or so, at least for
+ // non-turbocharged engines.
+ _compression = 8;
+ _displacement = power * (2*CIN2CM/HP2W);
}
void PistonEngine::setTurboParams(float turbo, float maxMP)
// This changes the "sea level" manifold air density
float P0 = Atmosphere::getStdPressure(0);
- float P = P0 * _turbo;
+ float P = P0 * (1 + _boost * (_turbo - 1));
if(P > _maxMP) P = _maxMP;
float T = Atmosphere::getStdTemperature(0) * Math::pow(P/P0, 2./7.);
_rho0 = P / (287.1 * T);
}
-float PistonEngine::getPower()
+void PistonEngine::setDisplacement(float d)
+{
+ _displacement = d;
+}
+
+void PistonEngine::setCompression(float c)
+{
+ _compression = c;
+}
+
+float PistonEngine::getMaxPower()
{
return _power0;
}
_mixture = m;
}
-void PistonEngine::calc(float P, float T, float speed,
- float* torqueOut, float* fuelFlowOut)
+void PistonEngine::setBoost(float boost)
+{
+ _boost = boost;
+}
+
+float PistonEngine::getTorque()
{
- // The actual fuel flow
- float fuel = _mixture * _mixCoeff * speed;
-
- // manifold air density
- if(_turbo != 1) {
- float P1 = P * _turbo;
- if(P1 > _maxMP) P1 = _maxMP;
- T *= Math::pow(P1/P, 2./7.);
- P = P1;
- }
- float density = P / (287.1 * T);
-
- float rho = density * _throttle;
+ return _torque;
+}
+
+float PistonEngine::getFuelFlow()
+{
+ return _fuelFlow;
+}
+
+float PistonEngine::getMP()
+{
+ return _mp;
+}
+
+float PistonEngine::getEGT()
+{
+ return _egt;
+}
+
+void PistonEngine::calc(float pressure, float temp, float speed)
+{
+ // Calculate manifold pressure as ambient pressure modified for
+ // turbocharging and reduced by the throttle setting. According
+ // to Dave Luff, minimum throttle at sea level corresponds to 6"
+ // manifold pressure. Assume that this means that minimum MP is
+ // always 20% of ambient pressure.
+ _mp = pressure * (1 + _boost*(_turbo-1)); // turbocharger
+ _mp *= (0.2 + 0.8 * _throttle); // throttle
+ if(_mp > _maxMP) _mp = _maxMP; // wastegate
+
+ // Air entering the manifold does so rapidly, and thus the
+ // pressure change can be assumed to be adiabatic. Calculate a
+ // temperature change, and use that to get the density.
+ float T = temp * Math::pow(_mp/pressure, 2.0/7.0);
+ float rho = _mp / (287.1 * T);
+
+ // The actual fuel flow is determined only by engine RPM and the
+ // mixture setting. Not all of this will burn with the same
+ // efficiency.
+ _fuelFlow = _mixture * speed * _mixCoeff;
// How much fuel could be burned with ideal (i.e. uncorrected!)
// combustion.
// interpolate. This vaguely matches a curve I copied out of a
// book for a single engine. Shrug.
float burned;
- float r = fuel/burnable;
+ float r = _fuelFlow/burnable;
if (burnable == 0) burned = 0;
- else if(r < .625) burned = fuel;
+ else if(r < .625) burned = _fuelFlow;
else if(r > 1.375) burned = burnable;
- else burned = fuel + (burnable-fuel)*(r-.625)*(4.0/3.0);
+ else
+ burned = _fuelFlow + (burnable-_fuelFlow)*(r-.625)*(4.0/3.0);
// And finally the power is just the reference power scaled by the
- // amount of fuel burned.
+ // amount of fuel burned, and torque is that divided by RPM.
float power = _power0 * burned/_f0;
-
- *torqueOut = power/speed;
- *fuelFlowOut = fuel;
+ _torque = power/speed;
+
+ // Now EGT. This one gets a little goofy. We can calculate the
+ // work done by an isentropically expanding exhaust gas as the
+ // mass of the gas times the specific heat times the change in
+ // temperature. The mass is just the engine displacement times
+ // the manifold density, plus the mass of the fuel, which we know.
+ // The change in temperature can be calculated adiabatically as a
+ // function of the exhaust gas temperature and the compression
+ // ratio (which we know). So just rearrange the equation to get
+ // EGT as a function of engine power. Cool. I'm using a value of
+ // 1300 J/(kg*K) for the exhaust gas specific heat. I found this
+ // on a web page somewhere; no idea if it's accurate. Also,
+ // remember that four stroke engines do one combustion cycle every
+ // TWO revolutions, so the displacement per revolution is half of
+ // what we'd expect. And diddle the work done by the gas a bit to
+ // account for non-thermodynamic losses like internal friction;
+ // 10% should do it.
+
+ float massFlow = _fuelFlow + (rho * 0.5 * _displacement * speed);
+ float specHeat = 1300;
+ float corr = 1.0/(Math::pow(_compression, 0.4) - 1);
+ _egt = corr * (power * 1.1) / (massFlow * specHeat);
}
}; // namespace yasim
// Initializes an engine from known "takeoff" parameters.
PistonEngine(float power, float spd);
void setTurboParams(float mul, float maxMP);
+ void setDisplacement(float d);
+ void setCompression(float c);
void setThrottle(float throttle);
void setMixture(float mixture);
+ void setBoost(float boost); // fraction of turbo-mul used
- float getPower();
+ float getMaxPower(); // max sea-level power
- // Calculates power output and fuel flow, based on a given
- // throttle setting (0-1 corresponding to the fraction of
- // "available" manifold pressure), mixture (fuel flux per rpm,
- // 0-1, where 1 is "max rich", or a little bit more than needed
- // for rated power at sea level)
- void calc(float pressure, float temp, float speed,
- float* powerOut, float* fuelFlowOut);
+ void calc(float pressure, float temp, float speed);
+ float getTorque();
+ float getFuelFlow();
+ float getMP();
+ float getEGT();
private:
+ // Static configuration:
float _power0; // reference power setting
float _omega0; // " engine speed
float _rho0; // " manifold air density
float _f0; // "ideal" fuel flow at P0/omega0
float _mixCoeff; // fuel flow per omega at full mixture
+ float _turbo; // (or super-)charger pressure multiplier
+ float _maxMP; // wastegate setting
+ float _displacement; // piston stroke volume
+ float _compression; // compression ratio (>1)
// Runtime settables:
float _throttle;
float _mixture;
+ float _boost;
- float _turbo;
- float _maxMP;
+ // Runtime state/output:
+ float _mp;
+ float _torque;
+ float _fuelFlow;
+ float _egt;
};
}; // namespace yasim
bool goingUp = false;
float step = 10;
while(true) {
- float etau, ptau, dummy;
+ float ptau, dummy;
_prop->calc(_rho, speed, _omega, &dummy, &ptau);
- _eng->calc(_pressure, _temp, _omega, &etau, &dummy);
+ _eng->calc(_pressure, _temp, _omega);
+ float etau = _eng->getTorque();
float tdiff = etau - ptau;
if(Math::abs(tdiff/_moment) < 0.1)
_eng->setThrottle(_throttle);
_eng->setMixture(_mixture);
- _prop->calc(_rho, speed, _omega,
- &thrust, &propTorque);
- _eng->calc(_pressure, _temp, _omega, &engTorque, &_fuelFlow);
+ _prop->calc(_rho, speed, _omega, &thrust, &propTorque);
+ _eng->calc(_pressure, _temp, _omega);
+ engTorque = _eng->getTorque();
+ _fuelFlow = _eng->getFuelFlow();
// Turn the thrust into a vector and save it
Math::mul3(thrust, _dir, _thrust);
--- /dev/null
+#include "Math.hpp"
+#include "SimpleJet.hpp"
+
+namespace yasim {
+
+SimpleJet::SimpleJet()
+{
+ _thrust = 0;
+}
+
+void SimpleJet::setThrust(float thrust)
+{
+ _thrust = thrust;
+}
+
+void SimpleJet::getThrust(float* out)
+{
+ Math::mul3(_thrust * _throttle, _dir, out);
+}
+
+void SimpleJet::getTorque(float* out)
+{
+ out[0] = out[1] = out[2] = 0;
+}
+
+void SimpleJet::getGyro(float* out)
+{
+ out[0] = out[1] = out[2] = 0;
+}
+
+float SimpleJet::getFuelFlow()
+{
+ return 0;
+}
+
+void SimpleJet::integrate(float dt)
+{
+ return;
+}
+
+void SimpleJet::stabilize()
+{
+ return;
+}
+
+}; // namespace yasim
--- /dev/null
+#ifndef _SIMPLEJET_HPP
+#define _SIMPLEJET_HPP
+
+#include "Thruster.hpp"
+
+namespace yasim {
+
+// As simple a Thruster subclass as you can find. It makes thrust. Period.
+class SimpleJet : public Thruster
+{
+public:
+ SimpleJet();
+ void setThrust(float thrust);
+ virtual void getThrust(float* out);
+ virtual void getTorque(float* out);
+ virtual void getGyro(float* out);
+ virtual float getFuelFlow();
+ virtual void integrate(float dt);
+ virtual void stabilize();
+private:
+ float _thrust;
+};
+
+}; // namespace yasim
+#endif // _SIMPLEJET_HPP
if(debugCount >= 3) {
debugCount = 0;
-// printf("RPM %.1f FF %.1f\n",
-// fgGetFloat("/engines/engine[0]/rpm"),
-// fgGetFloat("/engines/engine[0]/fuel-flow-gph"));
+// printf("N1 %5.1f N2 %5.1f FF %7.1f EPR %4.2f EGT %6.1f\n",
+// fgGetFloat("/engines/engine[0]/n1"),
+// fgGetFloat("/engines/engine[0]/n2"),
+// fgGetFloat("/engines/engine[0]/fuel-flow-gph"),
+// fgGetFloat("/engines/engine[0]/epr"),
+// fgGetFloat("/engines/engine[0]/egt"));
// printf("gear: %f\n", fgGetFloat("/controls/gear-down"));
set_delta_t(dt);
_fdm = new FGFDM();
- // Because the integration method is via fourth-order Runge-Kutta,
- // we only get an "output" state for every 4 times the internal
- // forces are calculated. So divide dt by four to account for
- // this, and only run an iteration every fourth time through
- // update.
- _dt = dt * 4;
+ _dt = dt;
+
_fdm->getAirplane()->getModel()->getIntegrator()->setInterval(_dt);
- _updateCount = 0;
}
void YASim::report()
}
}
+void YASim::bind()
+{
+ // Run the superclass bind to set up a bunch of property ties
+ FGInterface::bind();
+
+ // Now UNtie the ones that we are going to set ourselves.
+ fgUntie("/consumables/fuel/tank[0]/level-gal_us");
+ fgUntie("/consumables/fuel/tank[1]/level-gal_us");
+
+ char buf[256];
+ for(int i=0; i<_fdm->getAirplane()->getModel()->numThrusters(); i++) {
+ sprintf(buf, "/engines/engine[%d]/fuel-flow-gph", i); fgUntie(buf);
+ sprintf(buf, "/engines/engine[%d]/rpm", i); fgUntie(buf);
+ sprintf(buf, "/engines/engine[%d]/mp-osi", i); fgUntie(buf);
+ sprintf(buf, "/engines/engine[%d]/egt-degf", i); fgUntie(buf);
+ }
+
+}
+
void YASim::init()
{
Airplane* a = _fdm->getAirplane();
sprintf(buf, "/controls/propeller-pitch[%d]", i); fgSetFloat(buf, 1);
sprintf(buf, "/controls/afterburner[%d]", i); fgSetFloat(buf, 0);
}
-
+
+ fgSetFloat("/controls/slats", 0);
+ fgSetFloat("/controls/spoilers", 0);
// Are we at ground level? If so, lift the plane up so the gear
// clear the ground
// Blank the state, and copy in ours
State s;
m->setState(&s);
-
copyToYASim(true);
+
+ _fdm->getExternalInput();
+ _fdm->getAirplane()->stabilizeThrust();
+
set_inited(true);
}
int i;
for(i=0; i<iterations; i++) {
- // Remember, update only every 4th call
- _updateCount++;
- if(_updateCount >= 4) {
-
copyToYASim(false);
_fdm->iterate(_dt);
copyFromYASim();
printDEBUG();
-
- _updateCount = 0;
- }
}
}
pe->getTorque(tmp);
float power = Math::mag3(tmp) * pe->getOmega();
- float maxPower = pe->getPistonEngine()->getPower();
+ float maxPower = pe->getPistonEngine()->getMaxPower();
fge->set_MaxHP(maxPower * W2HP);
fge->set_Percentage_Power(100 * power/maxPower);
// Load externally set stuff into the FDM
virtual void init();
+ virtual void bind();
// Run an iteration
virtual void update(int iterations);
void printDEBUG();
yasim::FGFDM* _fdm;
- int _updateCount;
float _dt;
};