src/FDM/YASim/PistonEngine.cpp

   1 #include "Atmosphere.hpp"
   2 #include "Math.hpp"
   3 #include "PistonEngine.hpp"
   4 namespace yasim {
   5
   6 const static float HP2W = 745.7f;
   7 const static float CIN2CM = 1.6387064e-5f;
   8 const static float RPM2RADPS = 0.1047198f;
   9
  10 PistonEngine::PistonEngine(float power, float speed)
  11 {
  12     _boost = 1;
  13     _running = false;
  14     _cranking = false;
  15     _fuel = true;
  16
  17     // Presume a BSFC (in lb/hour per HP) of 0.45.  In SI that becomes
  18     // (2.2 lb/kg, 745.7 W/hp, 3600 sec/hour) 7.62e-08 kg/Ws.
  19     _f0 = power * 7.62e-08f;
  20
  21     _power0 = power;
  22     _omega0 = speed;
  23
  24     // We must be at sea level under standard conditions
  25     _rho0 = Atmosphere::getStdDensity(0);
  26
  27     // Further presume that takeoff is (duh) full throttle and
  28     // peak-power, that means that by our efficiency function, we are
  29     // at 11/8 of "ideal" fuel flow.
  30     float realFlow = _f0 * (11.0f/8.0f);
  31     _mixCoeff = realFlow * 1.1f / _omega0;
  32
  33     _turbo = 1;
  34     _maxMP = 1e6; // No waste gate on non-turbo engines.
  35
  36     // Guess at reasonable values for these guys.  Displacements run
  37     // at about 2 cubic inches per horsepower or so, at least for
  38     // non-turbocharged engines.
  39     _compression = 8;
  40     _displacement = power * (2*CIN2CM/HP2W);
  41 }
  42
  43 void PistonEngine::setTurboParams(float turbo, float maxMP)
  44 {
  45     _turbo = turbo;
  46     _maxMP = maxMP;
  47
  48     // This changes the "sea level" manifold air density
  49     float P0 = Atmosphere::getStdPressure(0);
  50     float P = P0 * (1 + _boost * (_turbo - 1));
  51     if(P > _maxMP) P = _maxMP;
  52     float T = Atmosphere::getStdTemperature(0) * Math::pow(P/P0, 2./7.);
  53     _rho0 = P / (287.1f * T);
  54 }
  55
  56 void PistonEngine::setDisplacement(float d)
  57 {
  58     _displacement = d;
  59 }
  60
  61 void PistonEngine::setCompression(float c)
  62 {
  63     _compression = c;
  64 }
  65
  66 float PistonEngine::getMaxPower()
  67 {
  68     return _power0;
  69 }
  70
  71 bool PistonEngine::isCranking()
  72 {
  73     return _cranking;
  74 }
  75
  76 float PistonEngine::getTorque()
  77 {
  78     return _torque;
  79 }
  80
  81 float PistonEngine::getFuelFlow()
  82 {
  83     return _fuelFlow;
  84 }
  85
  86 float PistonEngine::getMP()
  87 {
  88     return _mp;
  89 }
  90
  91 float PistonEngine::getEGT()
  92 {
  93     return _egt;
  94 }
  95
  96 void PistonEngine::calc(float pressure, float temp, float speed)
  97 {
  98     if(_magnetos == 0 || speed < 60*RPM2RADPS)
  99         _running = false;
 100     else if(_fuel == false)
 101         _running = false;
 102     else
 103         _running = true;
 104
 105     // Calculate manifold pressure as ambient pressure modified for
 106     // turbocharging and reduced by the throttle setting.  According
 107     // to Dave Luff, minimum throttle at sea level corresponds to 6"
 108     // manifold pressure.  Assume that this means that minimum MP is
 109     // always 20% of ambient pressure. (But that's too much idle
 110     // power, so use 10% instead!) But we need to produce _zero_
 111     // thrust at that setting, so hold onto the "output" value
 112     // separately.  Ick.
 113     _mp = pressure * (1 + _boost*(_turbo-1)); // turbocharger
 114     float mp = _mp * (0.1f + 0.9f * _throttle); // throttle
 115     _mp *= _throttle;
 116     if(mp > _maxMP) mp = _maxMP;              // wastegate
 117
 118     // Air entering the manifold does so rapidly, and thus the
 119     // pressure change can be assumed to be adiabatic.  Calculate a
 120     // temperature change, and use that to get the density.
 121     float T = temp * Math::pow(mp/pressure, 2.0/7.0);
 122     float rho = mp / (287.1f * T);
 123
 124     // The actual fuel flow is determined only by engine RPM and the
 125     // mixture setting.  Not all of this will burn with the same
 126     // efficiency.
 127     _fuelFlow = _mixture * speed * _mixCoeff;
 128     if(_fuel == false) _fuelFlow = 0;
 129
 130     // How much fuel could be burned with ideal (i.e. uncorrected!)
 131     // combustion.
 132     float burnable = _f0 * (rho/_rho0) * (speed/_omega0);
 133
 134     // Calculate the fuel that actually burns to produce work.  The
 135     // idea is that less than 5/8 of ideal, we get complete
 136     // combustion.  We use up all the oxygen at 1 3/8 of ideal (that
 137     // is, you need to waste fuel to use all your O2).  In between,
 138     // interpolate.  This vaguely matches a curve I copied out of a
 139     // book for a single engine.  Shrug.
 140     float burned;
 141     float r = _fuelFlow/burnable;
 142     if     (burnable == 0) burned = 0;
 143     else if(r < .625)      burned = _fuelFlow;
 144     else if(r > 1.375)     burned = burnable;
 145     else
 146         burned = _fuelFlow + (burnable-_fuelFlow)*(r-0.625f)*(4.0f/3.0f);
 147
 148     // Correct for engine control state
 149     if(!_running)
 150         burned = 0;
 151     if(_magnetos < 3)
 152         burned *= 0.9f;
 153
 154     // And finally the power is just the reference power scaled by the
 155     // amount of fuel burned, and torque is that divided by RPM.
 156     float power = _power0 * burned/_f0;
 157     _torque = power/speed;
 158
 159     // Figure that the starter motor produces 15% of the engine's
 160     // cruise torque.  Assuming 60RPM starter speed vs. 1800RPM cruise
 161     // speed on a 160HP engine, that comes out to about 160*.15/30 ==
 162     // 0.8 HP starter motor.  Which sounds about right to me.  I think
 163     // I've finally got this tuned. :)
 164     if(_cranking && !_running)
 165         _torque += 0.15f * _power0/_omega0;
 166
 167     // Also, add a negative torque of 8% of cruise, to represent
 168     // internal friction.  Propeller aerodynamic friction is too low
 169     // at low RPMs to provide a good deceleration.  Interpolate it
 170     // away as we approach cruise RPMs (full at 50%, zero at 100%),
 171     // though, to prevent interaction with the power computations.
 172     // Ugly.
 173     if(speed > 0 && speed < _omega0) {
 174         float interp = 2 - 2*speed/_omega0;
 175         interp = (interp > 1) ? 1 : interp;
 176         _torque -= 0.08f * (_power0/_omega0) * interp;
 177     }
 178
 179     // Now EGT.  This one gets a little goofy.  We can calculate the
 180     // work done by an isentropically expanding exhaust gas as the
 181     // mass of the gas times the specific heat times the change in
 182     // temperature.  The mass is just the engine displacement times
 183     // the manifold density, plus the mass of the fuel, which we know.
 184     // The change in temperature can be calculated adiabatically as a
 185     // function of the exhaust gas temperature and the compression
 186     // ratio (which we know).  So just rearrange the equation to get
 187     // EGT as a function of engine power.  Cool.  I'm using a value of
 188     // 1300 J/(kg*K) for the exhaust gas specific heat.  I found this
 189     // on a web page somewhere; no idea if it's accurate.  Also,
 190     // remember that four stroke engines do one combustion cycle every
 191     // TWO revolutions, so the displacement per revolution is half of
 192     // what we'd expect.  And diddle the work done by the gas a bit to
 193     // account for non-thermodynamic losses like internal friction;
 194     // 10% should do it.
 195
 196     float massFlow = _fuelFlow + (rho * 0.5f * _displacement * speed);
 197     float specHeat = 1300;
 198     float corr = 1.0f/(Math::pow(_compression, 0.4f) - 1.0f);
 199     _egt = corr * (power * 1.1f) / (massFlow * specHeat);
 200     if(_egt < temp) _egt = temp;
 201 }
 202
 203 }; // namespace yasim