+#include <stdio.h>
+
#include "Atmosphere.hpp"
#include "Math.hpp"
#include "Propeller.hpp"
namespace yasim {
Propeller::Propeller(float radius, float v, float omega,
- float rho, float power, float omega0,
- float power0)
+ float rho, float power)
{
// Initialize numeric constants:
- _lambdaPeak = Math::pow(9.0, -1.0/8.0);
- _beta = 1.0/(Math::pow(9.0, -1.0/8.0) - Math::pow(9.0, -9.0/8.0));
+ _lambdaPeak = Math::pow(5.0, -1.0/4.0);
+ _beta = 1.0f/(Math::pow(5.0f, -1.0f/4.0f) - Math::pow(5.0f, -5.0f/4.0f));
_r = radius;
- _etaC = 0.85; // make this settable?
+ _etaC = 0.85f; // make this settable?
- _J0 = v/(omega*_lambdaPeak);
+ _j0 = v/(omega*_lambdaPeak);
+ _baseJ0 = _j0;
float V2 = v*v + (_r*omega)*(_r*omega);
- _F0 = 2*_etaC*power/(rho*v*V2);
-
- float stallAngle = 0.25;
- _lambdaS = _r*(_J0/_r - stallAngle) / _J0;
-
- // Now compute a correction for zero forward speed to make the
- // takeoff performance correct.
- float torque0 = power0/omega0;
- float thrust, torque;
- _takeoffCoef = 1;
- calc(Atmosphere::getStdDensity(0), 0, omega0, &thrust, &torque);
- _takeoffCoef = torque/torque0;
+ _f0 = 2*_etaC*power/(rho*v*V2);
+
+ _matchTakeoff = false;
+ _manual = false;
+ _proppitch = 0;
+ _propfeather = 0;
+}
+
+void Propeller::setTakeoff(float omega0, float power0)
+{
+ // Takeoff thrust coefficient at lambda==0
+ _matchTakeoff = true;
+ float V2 = _r*omega0 * _r*omega0;
+ float gamma = _etaC * _beta / _j0;
+ float torque = power0 / omega0;
+ float density = Atmosphere::getStdDensity(0);
+ _tc0 = (torque * gamma) / (0.5f * density * V2 * _f0);
+}
+
+void Propeller::setStops(float fine_stop, float coarse_stop)
+{
+ _fine_stop = fine_stop;
+ _coarse_stop = coarse_stop;
+}
+
+void Propeller::modPitch(float mod)
+{
+ _j0 *= mod;
+ if(_j0 < _fine_stop*_baseJ0) _j0 = _fine_stop*_baseJ0;
+ if(_j0 > _coarse_stop*_baseJ0) _j0 = _coarse_stop*_baseJ0;
+}
+
+void Propeller::setManualPitch()
+{
+ _manual = true;
+}
+
+void Propeller::setPropPitch(float proppitch)
+{
+ // makes only positive range of axis effective.
+ _proppitch = Math::clamp(proppitch, 0, 1);
+}
+
+void Propeller::setPropFeather(int state)
+{
+ // 0 = normal, 1 = feathered
+ _propfeather = (state != 0);
}
void Propeller::calc(float density, float v, float omega,
float* thrustOut, float* torqueOut)
{
+ // For manual pitch, exponentially modulate the J0 value between
+ // 0.25 and 4. A prop pitch of 0.5 results in no change from the
+ // base value.
+ // TODO: integrate with _fine_stop and _coarse_stop variables
+ if (_manual)
+ _j0 = _baseJ0 * Math::pow(2, 2 - 4*_proppitch);
+
float tipspd = _r*omega;
float V2 = v*v + tipspd*tipspd;
- // Clamp v (forward velocity) to zero, now that we've used it to
- // calculate V (propeller "speed")
+ // Sanify
if(v < 0) v = 0;
+ if(omega < 0.001) omega = 0.001;
- float J = v/omega;
- float lambda = J/_J0;
+ float J = v/omega; // Advance ratio
+ float lambda = J/_j0; // Unitless scalar advance ratio
- float torque = 0;
- if(lambda > 1) {
- lambda = 1.0/lambda;
- torque = (density*V2*_F0*_J0)/(8*_etaC*_beta*(1-_lambdaPeak));
- }
+ // There's an undefined point at lambda == 1.
+ if(lambda == 1.0f) lambda = 0.9999f;
- // There's an undefined point at 1. Just offset by a tiny bit to
- // fix (note: the discontinuity is at EXACTLY one, this is about
- // the only time in history you'll see me use == on a floating
- // point number!)
- if(lambda == 1) lambda = 0.9999;
-
- // Compute thrust, remembering to clamp lambda to the stall point
- float lambda2 = lambda < _lambdaS ? _lambdaS : lambda;
- float thrust = (0.5*density*V2*_F0/(1-_lambdaPeak))*(1-lambda2);
-
- // Calculate lambda^8
- float l8 = lambda*lambda; l8 = l8*l8; l8 = l8*l8;
-
- // thrust/torque ratio
- float gamma = (_etaC*_beta/_J0)*(1-l8);
-
- // Correct slow speeds to get the takeoff parameters correct
- if(lambda < _lambdaPeak) {
- // This will interpolate takeoffCoef along a descending from 1
- // at lambda==0 to 0 at the peak, fairing smoothly into the
- // flat peak.
- float frac = (lambda - _lambdaPeak)/_lambdaPeak;
- gamma *= 1 + (_takeoffCoef - 1)*frac*frac;
- }
+ float l4 = lambda*lambda; l4 = l4*l4; // lambda^4
+ float gamma = (_etaC*_beta/_j0)*(1-l4); // thrust/torque ratio
+
+ // Compute a thrust coefficient, with clamping at very low
+ // lambdas (fast propeller / slow aircraft).
+ float tc = (1 - lambda) / (1 - _lambdaPeak);
+ if(_matchTakeoff && tc > _tc0) tc = _tc0;
- if(torque > 0) {
- torque -= thrust/gamma;
- thrust = -thrust;
- } else {
- torque = thrust/gamma;
+ float thrust = 0.5f * density * V2 * _f0 * tc;
+ float torque = thrust/gamma;
+ if(lambda > 1) {
+ // This is the negative thrust / windmilling regime. Throw
+ // out the efficiency graph approach and instead simply
+ // extrapolate the existing linear thrust coefficient and a
+ // torque coefficient that crosses the axis at a preset
+ // windmilling speed. The tau0 value is an analytically
+ // calculated (i.e. don't mess with it) value for a torque
+ // coefficient at lamda==1.
+ float tau0 = (0.25f * _j0) / (_etaC * _beta * (1 - _lambdaPeak));
+ float lambdaWM = 1.2f; // lambda of zero torque (windmilling)
+ torque = tau0 - tau0 * (lambda - 1) / (lambdaWM - 1);
+ torque *= 0.5f * density * V2 * _f0;
}
*thrustOut = thrust;