src/FDM/YASim/PropEngine.cpp

   1 #include "Math.hpp"
   2 #include "Propeller.hpp"
   3 #include "Engine.hpp"
   4 #include "PropEngine.hpp"
   5 namespace yasim {
   6
   7 PropEngine::PropEngine(Propeller* prop, Engine* eng, float moment)
   8 {
   9     // Start off at 500rpm, because the start code doesn't exist yet
  10     _omega = 52.3f;
  11     _dir[0] = 1; _dir[1] = 0; _dir[2] = 0;
  12
  13     _variable = false;
  14     _gearRatio = 1;
  15
  16     _prop = prop;
  17     _eng = eng;
  18     _moment = moment;
  19     _fuel = true;
  20     _contra = false;
  21 }
  22
  23 PropEngine::~PropEngine()
  24 {
  25     delete _prop;
  26     delete _eng;
  27 }
  28
  29 void PropEngine::setMagnetos(int pos)
  30 {
  31     _magnetos = pos;
  32 }
  33
  34 void PropEngine::setAdvance(float advance)
  35 {
  36     _advance = Math::clamp(advance, 0, 1);
  37 }
  38
  39 void PropEngine::setPropPitch(float proppitch)
  40 {
  41     // update Propeller property
  42     _prop->setPropPitch(proppitch);
  43 }
  44
  45 void PropEngine::setPropFeather(int state)
  46 {
  47     // toggle prop feathering on/off
  48     _prop->setPropFeather(state);
  49 }
  50
  51 void PropEngine::setVariableProp(float min, float max)
  52 {
  53     _variable = true;
  54     _minOmega = min;
  55     _maxOmega = max;
  56 }
  57
  58 bool PropEngine::isRunning()
  59 {
  60     return _eng->isRunning();
  61 }
  62
  63 bool PropEngine::isCranking()
  64 {
  65     return _eng->isCranking();
  66 }
  67
  68 float PropEngine::getOmega()
  69 {
  70     return _omega;
  71 }
  72
  73 void PropEngine::setOmega (float omega)
  74 {
  75     _omega = omega;
  76 }
  77
  78 void PropEngine::getThrust(float* out)
  79 {
  80     int i;
  81     for(i=0; i<3; i++) out[i] = _thrust[i];
  82 }
  83
  84 void PropEngine::getTorque(float* out)
  85 {
  86     int i;
  87     for(i=0; i<3; i++) out[i] = _torque[i];
  88 }
  89
  90 void PropEngine::getGyro(float* out)
  91 {
  92     int i;
  93     for(i=0; i<3; i++) out[i] = _gyro[i];
  94 }
  95
  96 float PropEngine::getFuelFlow()
  97 {
  98     return _fuelFlow;
  99 }
 100
 101 void PropEngine::stabilize()
 102 {
 103     float speed = -Math::dot3(_wind, _dir);
 104     _eng->setThrottle(_throttle);
 105     _eng->setMixture(_mixture);
 106
 107     _eng->setStarter(false);
 108     _eng->setMagnetos(3);
 109
 110     bool running_state = _eng->isRunning();
 111     _eng->setRunning(true);
 112
 113     if(_variable) {
 114         _omega = _minOmega + _advance * (_maxOmega - _minOmega);
 115         _prop->modPitch(1e6); // Start at maximum pitch and move down
 116     } else {
 117         _omega = 52;
 118     }
 119
 120     bool goingUp = false;
 121     float step = 10;
 122
 123     // If we cannot manage this in 100 iterations, give up.
 124     for (int n = 0; n < 100; n++) {
 125         float ptau, thrust;
 126         _prop->calc(_rho, speed, _omega * _gearRatio, &thrust, &ptau);
 127         _eng->calc(_pressure, _temp, _omega);
 128         _eng->stabilize();
 129
 130         // Do it again -- the turbo sets the target MP in the first
 131         // run, stabilize sets the current to the target, then we need
 132         // to run again to get the correct output torque.  Clumsy, but
 133         // it works without side effects (other than solver
 134         // performance).  In the future, the Engine objects should
 135         // store state to allow them to do the work themselves.
 136         _eng->calc(_pressure, _temp, _omega);
 137
 138         // Compute torque as seen by the engine's end of the gearbox.
 139         // The propeller will be moving more slowly (for gear ratios
 140         // less than one), so it's torque will be higher than the
 141         // engine's, so multiply by _gearRatio to get the engine-side
 142         // value.
 143         ptau *= _gearRatio;
 144         float etau = _eng->getTorque();
 145         float tdiff = etau - ptau;
 146
 147         Math::mul3(thrust, _dir, _thrust);
 148
 149         if(Math::abs(tdiff/(_moment * _gearRatio)) < 0.1)
 150             break;
 151
 152         if(tdiff > 0) {
 153             if(!goingUp) step *= 0.5f;
 154             goingUp = true;
 155             if(!_variable)  _omega += step;
 156             else            _prop->modPitch(1+(step*0.005f));
 157         } else {
 158             if(goingUp) step *= 0.5f;
 159             goingUp = false;
 160             if(!_variable)  _omega -= step;
 161             else            _prop->modPitch(1-(step*0.005f));
 162         }
 163     }
 164
 165     // ...and back off
 166     _eng->setRunning(running_state);
 167 }
 168
 169 void PropEngine::init()
 170 {
 171     _omega = 0.01f;
 172     _eng->setStarter(false);
 173     _eng->setMagnetos(0);
 174 }
 175
 176 void PropEngine::integrate(float dt)
 177 {
 178     float speed = -Math::dot3(_wind, _dir);
 179
 180     float propTorque, engTorque, thrust;
 181
 182     _eng->setThrottle(_throttle);
 183     _eng->setStarter(_starter);
 184     _eng->setMagnetos(_magnetos);
 185     _eng->setMixture(_mixture);
 186     _eng->setFuelState(_fuel);
 187
 188     _prop->calc(_rho, speed, _omega * _gearRatio, &thrust, &propTorque);
 189     if(_omega == 0.0)
 190         _omega = 0.001; // hack to get around reports of NaNs somewhere...
 191     propTorque *= _gearRatio;
 192     _eng->calc(_pressure, _temp, _omega);
 193     _eng->integrate(dt);
 194     engTorque = _eng->getTorque();
 195     _fuelFlow = _eng->getFuelFlow();
 196
 197     // Turn the thrust into a vector and save it
 198     Math::mul3(thrust, _dir, _thrust);
 199
 200     // We do our "RPM" computations on the engine's side of the
 201     // world, so modify the moment value accordingly.
 202     float momt = _moment * _gearRatio;
 203
 204     // Euler-integrate the RPM.  This doesn't need the full-on
 205     // Runge-Kutta stuff.
 206     float rotacc = (engTorque-propTorque)/Math::abs(momt);
 207     _omega += dt * rotacc;
 208     if (_omega < 0)
 209         _omega = 0 - _omega;    // don't allow negative RPM
 210                                 // FIXME: introduce proper windmilling
 211
 212     // Store the total angular momentum into _gyro, unless the
 213     // propeller is a counter-rotating pair (which has zero net
 214     // angular momentum, even though it *does* have an MoI for
 215     // acceleration purposes).
 216     Math::mul3(_contra ? 0 : _omega*momt, _dir, _gyro);
 217
 218     // Accumulate the engine torque, it acts on the body as a whole.
 219     // (Note: engine torque, not propeller torque.  They can be
 220     // different, but the difference goes to accelerating the
 221     // rotation.  It is the engine torque that is felt at the shaft
 222     // and works on the body.) (Note 2: contra-rotating propellers do
 223     // not exert net torque on the aircraft).
 224     float tau = _moment < 0 ? engTorque : -engTorque;
 225     Math::mul3(_contra ? 0 : tau, _dir, _torque);
 226
 227     // Iterate the propeller governor, if we have one.  Since engine
 228     // torque is basically constant with RPM, we want to make the
 229     // propeller torque at the target RPM equal to the engine by
 230     // varying the pitch.  Assume the the torque goes as the square of
 231     // the RPM (roughly correct) and compute a "target" torque for the
 232     // _current_ RPM.  Seek to that.  This is sort of a continuous
 233     // Newton-Raphson, basically.
 234     if(_variable) {
 235         float targetPropSpd = _minOmega + _advance*(_maxOmega-_minOmega);
 236         float targetOmega = targetPropSpd / _gearRatio; // -> "engine omega"
 237         float ratio2 = (_omega*_omega)/(targetOmega*targetOmega);
 238         float targetTorque = engTorque * ratio2;
 239
 240         float mod = propTorque < targetTorque ? 1.04f : (1.0f/1.04f);
 241
 242         // Convert to an acceleration here, so that big propellers
 243         // don't seek faster than small ones.
 244         float diff = Math::abs((propTorque - targetTorque) / momt);
 245         if(diff < 10) mod = 1 + (mod-1)*(0.1f*diff);
 246
 247         _prop->modPitch(mod);
 248     }
 249 }
 250
 251 }; // namespace yasim