src/FDM/YASim/Propeller.cpp

   1 #include "Atmosphere.hpp"
   2 #include "Math.hpp"
   3 #include "Propeller.hpp"
   4 namespace yasim {
   5
   6 Propeller::Propeller(float radius, float v, float omega,
   7                      float rho, float power, float omega0,
   8                      float power0)
   9 {
  10     // Initialize numeric constants:
  11     _lambdaPeak = Math::pow(9.0, -1.0/8.0);
  12     _beta = 1.0/(Math::pow(9.0, -1.0/8.0) - Math::pow(9.0, -9.0/8.0));
  13
  14     _r = radius;
  15     _etaC = 0.85; // make this settable?
  16
  17     _J0 = v/(omega*_lambdaPeak);
  18
  19     float V2 = v*v + (_r*omega)*(_r*omega);
  20     _F0 = 2*_etaC*power/(rho*v*V2);
  21
  22     float stallAngle = 0.25;
  23     _lambdaS = _r*(_J0/_r - stallAngle) / _J0;
  24
  25     // Now compute a correction for zero forward speed to make the
  26     // takeoff performance correct.
  27     float torque0 = power0/omega0;
  28     float thrust, torque;
  29     _takeoffCoef = 1;
  30     calc(Atmosphere::getStdDensity(0), 0, omega0, &thrust, &torque);
  31     _takeoffCoef = torque/torque0;
  32 }
  33
  34 void Propeller::calc(float density, float v, float omega,
  35                      float* thrustOut, float* torqueOut)
  36 {
  37     float tipspd = _r*omega;
  38     float V2 = v*v + tipspd*tipspd;
  39
  40     // Clamp v (forward velocity) to zero, now that we've used it to
  41     // calculate V (propeller "speed")
  42     if(v < 0) v = 0;
  43
  44     float J = v/omega;
  45     float lambda = J/_J0;
  46
  47     float torque = 0;
  48     if(lambda > 1) {
  49         lambda = 1.0/lambda;
  50         torque = (density*V2*_F0*_J0)/(8*_etaC*_beta*(1-_lambdaPeak));
  51     }
  52
  53     // There's an undefined point at 1.  Just offset by a tiny bit to
  54     // fix (note: the discontinuity is at EXACTLY one, this is about
  55     // the only time in history you'll see me use == on a floating
  56     // point number!)
  57     if(lambda == 1) lambda = 0.9999;
  58
  59     // Compute thrust, remembering to clamp lambda to the stall point
  60     float lambda2 = lambda < _lambdaS ? _lambdaS : lambda;
  61     float thrust = (0.5*density*V2*_F0/(1-_lambdaPeak))*(1-lambda2);
  62
  63     // Calculate lambda^8
  64     float l8 = lambda*lambda; l8 = l8*l8; l8 = l8*l8;
  65
  66     // thrust/torque ratio
  67     float gamma = (_etaC*_beta/_J0)*(1-l8);
  68
  69     // Correct slow speeds to get the takeoff parameters correct
  70     if(lambda < _lambdaPeak) {
  71         // This will interpolate takeoffCoef along a descending from 1
  72         // at lambda==0 to 0 at the peak, fairing smoothly into the
  73         // flat peak.
  74         float frac = (lambda - _lambdaPeak)/_lambdaPeak;
  75         gamma *= 1 + (_takeoffCoef - 1)*frac*frac;
  76     }
  77
  78     if(torque > 0) {
  79         torque -= thrust/gamma;
  80         thrust = -thrust;
  81     } else {
  82         torque = thrust/gamma;
  83     }
  84
  85     *thrustOut = thrust;
  86     *torqueOut = torque;
  87 }
  88
  89 }; // namespace yasim