src/FDM/JSBSim/models/propulsion/FGPiston.cpp

   1 /*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   2
   3  Module:       FGPiston.cpp
   4  Author:       Jon S. Berndt, JSBSim framework
   5                Dave Luff, Piston engine model
   6                Ronald Jensen, Piston engine model
   7  Date started: 09/12/2000
   8  Purpose:      This module models a Piston engine
   9
  10  ------------- Copyright (C) 2000  Jon S. Berndt (jsb@hal-pc.org) --------------
  11
  12  This program is free software; you can redistribute it and/or modify it under
  13  the terms of the GNU Lesser General Public License as published by the Free Software
  14  Foundation; either version 2 of the License, or (at your option) any later
  15  version.
  16
  17  This program is distributed in the hope that it will be useful, but WITHOUT
  18  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
  19  FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public License for more
  20  details.
  21
  22  You should have received a copy of the GNU Lesser General Public License along with
  23  this program; if not, write to the Free Software Foundation, Inc., 59 Temple
  24  Place - Suite 330, Boston, MA  02111-1307, USA.
  25
  26  Further information about the GNU Lesser General Public License can also be found on
  27  the world wide web at http://www.gnu.org.
  28
  29 FUNCTIONAL DESCRIPTION
  30 --------------------------------------------------------------------------------
  31
  32 This class descends from the FGEngine class and models a Piston engine based on
  33 parameters given in the engine config file for this class
  34
  35 HISTORY
  36 --------------------------------------------------------------------------------
  37 09/12/2000  JSB  Created
  38
  39 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  40 INCLUDES
  41 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
  42
  43 #include <sstream>
  44
  45 #include "FGPiston.h"
  46 #include <models/FGPropulsion.h>
  47 #include "FGPropeller.h"
  48
  49 namespace JSBSim {
  50
  51 static const char *IdSrc = "$Id$";
  52 static const char *IdHdr = ID_PISTON;
  53
  54 /*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  55 CLASS IMPLEMENTATION
  56 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
  57
  58 FGPiston::FGPiston(FGFDMExec* exec, Element* el, int engine_number)
  59   : FGEngine(exec, el, engine_number),
  60   R_air(287.3),                  // Gas constant for air J/Kg/K
  61   rho_fuel(800),                 // estimate
  62   calorific_value_fuel(47.3e6),
  63   Cp_air(1005),                  // Specific heat (constant pressure) J/Kg/K
  64   Cp_fuel(1700)
  65 {
  66   string token;
  67
  68   // Defaults and initializations
  69
  70   Type = etPiston;
  71   dt = State->Getdt();
  72
  73   // These items are read from the configuration file
  74
  75   Cycles = 2;
  76   IdleRPM = 600;
  77   MaxRPM = 2800;
  78   Displacement = 360;
  79   SparkFailDrop = 1.0;
  80   MaxHP = 200;
  81   MinManifoldPressure_inHg = 6.5;
  82   MaxManifoldPressure_inHg = 28.5;
  83   BSFC = 0.45;
  84
  85   // These are internal program variables
  86
  87   crank_counter = 0;
  88   Magnetos = 0;
  89   minMAP = 21950;
  90   maxMAP = 96250;
  91
  92   ResetToIC();
  93
  94   // Supercharging
  95   BoostSpeeds = 0;  // Default to no supercharging
  96   BoostSpeed = 0;
  97   Boosted = false;
  98   BoostOverride = 0;
  99   bBoostOverride = false;
 100   bTakeoffBoost = false;
 101   TakeoffBoost = 0.0;   // Default to no extra takeoff-boost
 102   int i;
 103   for (i=0; i<FG_MAX_BOOST_SPEEDS; i++) {
 104     RatedBoost[i] = 0.0;
 105     RatedPower[i] = 0.0;
 106     RatedAltitude[i] = 0.0;
 107     BoostMul[i] = 1.0;
 108     RatedMAP[i] = 100000;
 109     RatedRPM[i] = 2500;
 110     TakeoffMAP[i] = 100000;
 111   }
 112   for (i=0; i<FG_MAX_BOOST_SPEEDS-1; i++) {
 113     BoostSwitchAltitude[i] = 0.0;
 114     BoostSwitchPressure[i] = 0.0;
 115   }
 116   // Initialisation
 117   volumetric_efficiency = 0.8;  // Actually f(speed, load) but this will get us running
 118
 119   // First column is thi, second is neta (combustion efficiency)
 120   Lookup_Combustion_Efficiency = new FGTable(12);
 121   *Lookup_Combustion_Efficiency << 0.00 << 0.980;
 122   *Lookup_Combustion_Efficiency << 0.90 << 0.980;
 123   *Lookup_Combustion_Efficiency << 1.00 << 0.970;
 124   *Lookup_Combustion_Efficiency << 1.05 << 0.950;
 125   *Lookup_Combustion_Efficiency << 1.10 << 0.900;
 126   *Lookup_Combustion_Efficiency << 1.15 << 0.850;
 127   *Lookup_Combustion_Efficiency << 1.20 << 0.790;
 128   *Lookup_Combustion_Efficiency << 1.30 << 0.700;
 129   *Lookup_Combustion_Efficiency << 1.40 << 0.630;
 130   *Lookup_Combustion_Efficiency << 1.50 << 0.570;
 131   *Lookup_Combustion_Efficiency << 1.60 << 0.525;
 132   *Lookup_Combustion_Efficiency << 2.00 << 0.345;
 133
 134   Power_Mixture_Correlation = new FGTable(13);
 135   *Power_Mixture_Correlation << (14.7/1.6) << 0.780;
 136   *Power_Mixture_Correlation << 10 <<  0.860;
 137   *Power_Mixture_Correlation << 11 <<  0.935;
 138   *Power_Mixture_Correlation << 12 <<  0.980;
 139   *Power_Mixture_Correlation << 13 <<  1.000;
 140   *Power_Mixture_Correlation << 14 <<  0.990;
 141   *Power_Mixture_Correlation << 15 <<  0.964;
 142   *Power_Mixture_Correlation << 16 <<  0.925;
 143   *Power_Mixture_Correlation << 17 <<  0.880;
 144   *Power_Mixture_Correlation << 18 <<  0.830;
 145   *Power_Mixture_Correlation << 19 <<  0.785;
 146   *Power_Mixture_Correlation << 20 <<  0.740;
 147   *Power_Mixture_Correlation << (14.7/0.6) << 0.58;
 148
 149   Mixture_Efficiency_Correlation = new FGTable(15);
 150   *Mixture_Efficiency_Correlation << 0.05000 << 0.00000;
 151   *Mixture_Efficiency_Correlation << 0.05137 << 0.00862;
 152   *Mixture_Efficiency_Correlation << 0.05179 << 0.21552;
 153   *Mixture_Efficiency_Correlation << 0.05430 << 0.48276;
 154   *Mixture_Efficiency_Correlation << 0.05842 << 0.70690;
 155   *Mixture_Efficiency_Correlation << 0.06312 << 0.83621;
 156   *Mixture_Efficiency_Correlation << 0.06942 << 0.93103;
 157   *Mixture_Efficiency_Correlation << 0.07786 << 1.00000;
 158   *Mixture_Efficiency_Correlation << 0.08845 << 1.00000;
 159   *Mixture_Efficiency_Correlation << 0.09270 << 0.98276;
 160   *Mixture_Efficiency_Correlation << 0.10120 << 0.93103;
 161   *Mixture_Efficiency_Correlation << 0.11455 << 0.72414;
 162   *Mixture_Efficiency_Correlation << 0.12158 << 0.45690;
 163   *Mixture_Efficiency_Correlation << 0.12435 << 0.23276;
 164   *Mixture_Efficiency_Correlation << 0.12500 << 0.00000;
 165
 166
 167 /*
 168 Manifold_Pressure_Lookup = new
 169
 170         0       0.2      0.4      0.6      0.8     1
 171 0        1.0000    1.0000    1.0000    1.0000    1.0000    1.0000
 172 1000    0.7778    0.8212    0.8647    0.9081    0.9516    0.9950
 173 2000    0.5556    0.6424    0.7293    0.8162    0.9031    0.9900
 174 3000    0.3333    0.4637    0.5940    0.7243    0.8547    0.9850
 175 4000    0.2000    0.2849    0.4587    0.6324    0.8062    0.9800
 176 5000    0.2000    0.2000    0.3233    0.5406    0.7578    0.9750
 177 6000    0.2000    0.2000    0.2000    0.4487    0.7093    0.9700
 178 7000    0.2000    0.2000    0.2000    0.2000    0.4570    0.7611
 179 8000    0.2000    0.2000    0.2000    0.2000    0.2047    0.5522
 180 */
 181
 182   // Read inputs from engine data file where present.
 183
 184   if (el->FindElement("minmp")) // Should have ELSE statement telling default value used?
 185     MinManifoldPressure_inHg = el->FindElementValueAsNumberConvertTo("minmp","INHG");
 186   if (el->FindElement("maxmp"))
 187     MaxManifoldPressure_inHg = el->FindElementValueAsNumberConvertTo("maxmp","INHG");
 188   if (el->FindElement("displacement"))
 189     Displacement = el->FindElementValueAsNumberConvertTo("displacement","IN3");
 190   if (el->FindElement("maxhp"))
 191     MaxHP = el->FindElementValueAsNumberConvertTo("maxhp","HP");
 192   if (el->FindElement("sparkfaildrop"))
 193     SparkFailDrop = Constrain(0, 1 - el->FindElementValueAsNumber("sparkfaildrop"), 1);
 194   if (el->FindElement("cycles"))
 195     Cycles = el->FindElementValueAsNumber("cycles");
 196   if (el->FindElement("idlerpm"))
 197     IdleRPM = el->FindElementValueAsNumber("idlerpm");
 198   if (el->FindElement("maxrpm"))
 199     MaxRPM = el->FindElementValueAsNumber("maxrpm");
 200   if (el->FindElement("maxthrottle"))
 201     MaxThrottle = el->FindElementValueAsNumber("maxthrottle");
 202   if (el->FindElement("minthrottle"))
 203     MinThrottle = el->FindElementValueAsNumber("minthrottle");
 204   if (el->FindElement("bsfc"))
 205     BSFC = el->FindElementValueAsNumber("bsfc");
 206   if (el->FindElement("numboostspeeds")) { // Turbo- and super-charging parameters
 207     BoostSpeeds = (int)el->FindElementValueAsNumber("numboostspeeds");
 208     if (el->FindElement("boostoverride"))
 209       BoostOverride = (int)el->FindElementValueAsNumber("boostoverride");
 210     if (el->FindElement("takeoffboost"))
 211       TakeoffBoost = el->FindElementValueAsNumberConvertTo("takeoffboost", "PSI");
 212     if (el->FindElement("ratedboost1"))
 213       RatedBoost[0] = el->FindElementValueAsNumberConvertTo("ratedboost1", "PSI");
 214     if (el->FindElement("ratedboost2"))
 215       RatedBoost[1] = el->FindElementValueAsNumberConvertTo("ratedboost2", "PSI");
 216     if (el->FindElement("ratedboost3"))
 217       RatedBoost[2] = el->FindElementValueAsNumberConvertTo("ratedboost3", "PSI");
 218     if (el->FindElement("ratedpower1"))
 219       RatedPower[0] = el->FindElementValueAsNumberConvertTo("ratedpower1", "HP");
 220     if (el->FindElement("ratedpower2"))
 221       RatedPower[1] = el->FindElementValueAsNumberConvertTo("ratedpower2", "HP");
 222     if (el->FindElement("ratedpower3"))
 223       RatedPower[2] = el->FindElementValueAsNumberConvertTo("ratedpower3", "HP");
 224     if (el->FindElement("ratedrpm1"))
 225       RatedRPM[0] = el->FindElementValueAsNumber("ratedrpm1");
 226     if (el->FindElement("ratedrpm2"))
 227       RatedRPM[1] = el->FindElementValueAsNumber("ratedrpm2");
 228     if (el->FindElement("ratedrpm3"))
 229       RatedRPM[2] = el->FindElementValueAsNumber("ratedrpm3");
 230     if (el->FindElement("ratedaltitude1"))
 231       RatedAltitude[0] = el->FindElementValueAsNumberConvertTo("ratedaltitude1", "FT");
 232     if (el->FindElement("ratedaltitude2"))
 233       RatedAltitude[1] = el->FindElementValueAsNumberConvertTo("ratedaltitude2", "FT");
 234     if (el->FindElement("ratedaltitude3"))
 235       RatedAltitude[2] = el->FindElementValueAsNumberConvertTo("ratedaltitude3", "FT");
 236   }
 237   char property_name[80];
 238   snprintf(property_name, 80, "/engines/engine[%d]/power_hp", engine_number);
 239   PropertyManager->Tie(property_name, &HP);
 240   snprintf(property_name, 80, "/engines/engine[%d]/bsfc", engine_number);
 241   PropertyManager->Tie(property_name, &BSFC);
 242
 243   minMAP = MinManifoldPressure_inHg * inhgtopa;  // inHg to Pa
 244   maxMAP = MaxManifoldPressure_inHg * inhgtopa;
 245   StarterHP = sqrt(MaxHP) * 0.4;
 246
 247   // Set up and sanity-check the turbo/supercharging configuration based on the input values.
 248   if (TakeoffBoost > RatedBoost[0]) bTakeoffBoost = true;
 249   for (i=0; i<BoostSpeeds; ++i) {
 250     bool bad = false;
 251     if (RatedBoost[i] <= 0.0) bad = true;
 252     if (RatedPower[i] <= 0.0) bad = true;
 253     if (RatedAltitude[i] < 0.0) bad = true;  // 0.0 is deliberately allowed - this corresponds to unregulated supercharging.
 254     if (i > 0 && RatedAltitude[i] < RatedAltitude[i - 1]) bad = true;
 255     if (bad) {
 256       // We can't recover from the above - don't use this supercharger speed.
 257       BoostSpeeds--;
 258       // TODO - put out a massive error message!
 259       break;
 260     }
 261     // Now sanity-check stuff that is recoverable.
 262     if (i < BoostSpeeds - 1) {
 263       if (BoostSwitchAltitude[i] < RatedAltitude[i]) {
 264         // TODO - put out an error message
 265         // But we can also make a reasonable estimate, as below.
 266         BoostSwitchAltitude[i] = RatedAltitude[i] + 1000;
 267       }
 268       BoostSwitchPressure[i] = Atmosphere->GetPressure(BoostSwitchAltitude[i]) * psftopa;
 269       //cout << "BoostSwitchAlt = " << BoostSwitchAltitude[i] << ", pressure = " << BoostSwitchPressure[i] << '\n';
 270       // Assume there is some hysteresis on the supercharger gear switch, and guess the value for now
 271       BoostSwitchHysteresis = 1000;
 272     }
 273     // Now work out the supercharger pressure multiplier of this speed from the rated boost and altitude.
 274     RatedMAP[i] = Atmosphere->GetPressureSL() * psftopa + RatedBoost[i] * 6895;  // psi*6895 = Pa.
 275     // Sometimes a separate BCV setting for takeoff or extra power is fitted.
 276     if (TakeoffBoost > RatedBoost[0]) {
 277       // Assume that the effect on the BCV is the same whichever speed is in use.
 278       TakeoffMAP[i] = RatedMAP[i] + ((TakeoffBoost - RatedBoost[0]) * 6895);
 279       bTakeoffBoost = true;
 280     } else {
 281       TakeoffMAP[i] = RatedMAP[i];
 282       bTakeoffBoost = false;
 283     }
 284     BoostMul[i] = RatedMAP[i] / (Atmosphere->GetPressure(RatedAltitude[i]) * psftopa);
 285
 286   }
 287
 288   if (BoostSpeeds > 0) {
 289     Boosted = true;
 290     BoostSpeed = 0;
 291   }
 292   bBoostOverride = (BoostOverride == 1 ? true : false);
 293   if (MinThrottle < 0.001) MinThrottle = 0.001;  //MinThrottle is a denominator in a power equation so it can't be zero
 294   Debug(0); // Call Debug() routine from constructor if needed
 295 }
 296
 297 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 298
 299 FGPiston::~FGPiston()
 300 {
 301   char property_name[80];
 302   snprintf(property_name, 80, "/engines/engine[%d]/power_hp", EngineNumber);
 303   PropertyManager->Untie(property_name);
 304   snprintf(property_name, 80, "/engines/engine[%d]/bsfc", EngineNumber);
 305   PropertyManager->Untie(property_name);
 306
 307   delete Lookup_Combustion_Efficiency;
 308   delete Power_Mixture_Correlation;
 309   delete Mixture_Efficiency_Correlation;
 310   Debug(1); // Call Debug() routine from constructor if needed
 311 }
 312
 313 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 314
 315 void FGPiston::ResetToIC(void)
 316 {
 317   FGEngine::ResetToIC();
 318
 319   ManifoldPressure_inHg = Atmosphere->GetPressure() * psftoinhg; // psf to in Hg
 320   MAP = Atmosphere->GetPressure() * psftopa;
 321   double airTemperature_degK = RankineToKelvin(Atmosphere->GetTemperature());
 322   OilTemp_degK = airTemperature_degK;
 323   CylinderHeadTemp_degK = airTemperature_degK;
 324   ExhaustGasTemp_degK = airTemperature_degK;
 325   EGT_degC = ExhaustGasTemp_degK - 273;
 326   Thruster->SetRPM(0.0);
 327   RPM = 0.0;
 328 }
 329
 330 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 331
 332 double FGPiston::Calculate(void)
 333 {
 334   if (FuelFlow_gph > 0.0) ConsumeFuel();
 335
 336   Throttle = FCS->GetThrottlePos(EngineNumber);
 337   ThrottlePos = MinThrottle+((MaxThrottle-MinThrottle)*Throttle );
 338   Mixture = FCS->GetMixturePos(EngineNumber);
 339
 340   //
 341   // Input values.
 342   //
 343
 344   p_amb = Atmosphere->GetPressure() * psftopa;
 345   p_amb_sea_level = Atmosphere->GetPressureSL() * psftopa;
 346   T_amb = RankineToKelvin(Atmosphere->GetTemperature());
 347
 348   RPM = Thruster->GetRPM() * Thruster->GetGearRatio();
 349
 350   IAS = Auxiliary->GetVcalibratedKTS();
 351
 352   doEngineStartup();
 353   if (Boosted) doBoostControl();
 354   doMAP();
 355   doAirFlow();
 356   doFuelFlow();
 357
 358   //Now that the fuel flow is done check if the mixture is too lean to run the engine
 359   //Assume lean limit at 22 AFR for now - thats a thi of 0.668
 360   //This might be a bit generous, but since there's currently no audiable warning of impending
 361   //cutout in the form of misfiring and/or rough running its probably reasonable for now.
 362 //  if (equivalence_ratio < 0.668)
 363 //    Running = false;
 364
 365   doEnginePower();
 366 if(HP<0.1250)
 367   Running = false;
 368
 369   doEGT();
 370   doCHT();
 371   doOilTemperature();
 372   doOilPressure();
 373
 374   if (Thruster->GetType() == FGThruster::ttPropeller) {
 375     ((FGPropeller*)Thruster)->SetAdvance(FCS->GetPropAdvance(EngineNumber));
 376     ((FGPropeller*)Thruster)->SetFeather(FCS->GetPropFeather(EngineNumber));
 377   }
 378
 379   PowerAvailable = (HP * hptoftlbssec) - Thruster->GetPowerRequired();
 380
 381   return Thruster->Calculate(PowerAvailable);
 382 }
 383
 384 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 385
 386 double FGPiston::CalcFuelNeed(void)
 387 {
 388   return FuelFlow_gph / 3600 * 6 * State->Getdt() * Propulsion->GetRate();
 389 }
 390
 391 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 392
 393 int FGPiston::InitRunning(void) {
 394   Magnetos=3;
 395   //Thruster->SetRPM( 1.1*IdleRPM/Thruster->GetGearRatio() );
 396   Thruster->SetRPM( 1000 );
 397   Running=true;
 398   return 1;
 399 }
 400
 401 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 402 /**
 403  * Start or stop the engine.
 404  */
 405
 406 void FGPiston::doEngineStartup(void)
 407 {
 408   // Check parameters that may alter the operating state of the engine.
 409   // (spark, fuel, starter motor etc)
 410   bool spark;
 411   bool fuel;
 412
 413   // Check for spark
 414   Magneto_Left = false;
 415   Magneto_Right = false;
 416   // Magneto positions:
 417   // 0 -> off
 418   // 1 -> left only
 419   // 2 -> right only
 420   // 3 -> both
 421   if (Magnetos != 0) {
 422     spark = true;
 423   } else {
 424     spark = false;
 425   }  // neglects battery voltage, master on switch, etc for now.
 426
 427   if ((Magnetos == 1) || (Magnetos > 2)) Magneto_Left = true;
 428   if (Magnetos > 1)  Magneto_Right = true;
 429
 430   // Assume we have fuel for now
 431   fuel = !Starved;
 432
 433   // Check if we are turning the starter motor
 434   if (Cranking != Starter) {
 435     // This check saves .../cranking from getting updated every loop - they
 436     // only update when changed.
 437     Cranking = Starter;
 438     crank_counter = 0;
 439   }
 440
 441   if (Cranking) crank_counter++;  //Check mode of engine operation
 442
 443   if (!Running && spark && fuel) {  // start the engine if revs high enough
 444     if (Cranking) {
 445       if ((RPM > IdleRPM*0.8) && (crank_counter > 175)) // Add a little delay to startup
 446         Running = true;                         // on the starter
 447     } else {
 448       if (RPM > IdleRPM*0.8)                            // This allows us to in-air start
 449         Running = true;                         // when windmilling
 450     }
 451   }
 452
 453   // Cut the engine *power* - Note: the engine may continue to
 454   // spin if the prop is in a moving airstream
 455
 456   if ( Running && (!spark || !fuel) ) Running = false;
 457
 458   // Check for stalling (RPM = 0).
 459   if (Running) {
 460     if (RPM == 0) {
 461       Running = false;
 462     } else if ((RPM <= IdleRPM *0.8 ) && (Cranking)) {
 463       Running = false;
 464     }
 465   }
 466 }
 467
 468 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 469
 470 /**
 471  * Calculate the Current Boost Speed
 472  *
 473  * This function calculates the current turbo/supercharger boost speed
 474  * based on altitude and the (automatic) boost-speed control valve configuration.
 475  *
 476  * Inputs: p_amb, BoostSwitchPressure, BoostSwitchHysteresis
 477  *
 478  * Outputs: BoostSpeed
 479  */
 480
 481 void FGPiston::doBoostControl(void)
 482 {
 483   if(BoostSpeed < BoostSpeeds - 1) {
 484     // Check if we need to change to a higher boost speed
 485     if(p_amb < BoostSwitchPressure[BoostSpeed] - BoostSwitchHysteresis) {
 486       BoostSpeed++;
 487     }
 488   } else if(BoostSpeed > 0) {
 489     // Check if we need to change to a lower boost speed
 490     if(p_amb > BoostSwitchPressure[BoostSpeed - 1] + BoostSwitchHysteresis) {
 491       BoostSpeed--;
 492     }
 493   }
 494 }
 495
 496 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 497
 498 /**
 499  * Calculate the manifold absolute pressure (MAP) in inches hg
 500  *
 501  * This function calculates manifold absolute pressure (MAP)
 502  * from the throttle position, turbo/supercharger boost control
 503  * system, engine speed and local ambient air density.
 504  *
 505  * TODO: changes in MP should not be instantaneous -- introduce
 506  * a lag between throttle changes and MP changes, to allow pressure
 507  * to build up or disperse.
 508  *
 509  * Inputs: minMAP, maxMAP, p_amb, Throttle
 510  *
 511  * Outputs: MAP, ManifoldPressure_inHg
 512  */
 513
 514 void FGPiston::doMAP(void)
 515 {
 516     suction_loss = pow( ThrottlePos*0.98, RPM/MaxRPM );
 517     MAP = p_amb * suction_loss;
 518
 519     if(Boosted) {
 520       // If takeoff boost is fitted, we currently assume the following throttle map:
 521       // (In throttle % - actual input is 0 -> 1)
 522       // 99 / 100 - Takeoff boost
 523       // 96 / 97 / 98 - Rated boost
 524       // 0 - 95 - Idle to Rated boost (MinManifoldPressure to MaxManifoldPressure)
 525       // In real life, most planes would be fitted with a mechanical 'gate' between
 526       // the rated boost and takeoff boost positions.
 527       double T = Throttle; // processed throttle value.
 528       bool bTakeoffPos = false;
 529       if(bTakeoffBoost) {
 530         if(Throttle > 0.98) {
 531           //cout << "Takeoff Boost!!!!\n";
 532           bTakeoffPos = true;
 533         } else if(Throttle <= 0.95) {
 534           bTakeoffPos = false;
 535           T *= 1.0 / 0.95;
 536         } else {
 537           bTakeoffPos = false;
 538           //cout << "Rated Boost!!\n";
 539           T = 1.0;
 540         }
 541       }
 542       // Boost the manifold pressure.
 543       MAP += MAP * BoostMul[BoostSpeed] * RPM/MaxRPM;
 544       // Now clip the manifold pressure to BCV or Wastegate setting.
 545       if(bTakeoffPos) {
 546         if(MAP > TakeoffMAP[BoostSpeed]) {
 547           MAP = TakeoffMAP[BoostSpeed];
 548         }
 549       } else {
 550         if(MAP > RatedMAP[BoostSpeed]) {
 551           MAP = RatedMAP[BoostSpeed];
 552         }
 553       }
 554     }
 555
 556   // And set the value in American units as well
 557   ManifoldPressure_inHg = MAP / inhgtopa;
 558 }
 559
 560 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 561 /**
 562  * Calculate the air flow through the engine.
 563  * Also calculates ambient air density
 564  * (used in CHT calculation for air-cooled engines).
 565  *
 566  * Inputs: p_amb, R_air, T_amb, MAP, Displacement,
 567  *   RPM, volumetric_efficiency, ThrottlePos
 568  *
 569  * TODO: Model inlet manifold air temperature.
 570  *
 571  * Outputs: rho_air, m_dot_air
 572  */
 573
 574 void FGPiston::doAirFlow(void)
 575 {
 576
 577 rho_air = p_amb / (R_air * T_amb);
 578   double displacement_SI = Displacement * in3tom3;
 579   double swept_volume = (displacement_SI * (RPM/60)) / 2;
 580   double v_dot_air = swept_volume * volumetric_efficiency;
 581
 582   double rho_air_manifold = MAP / (R_air * T_amb);
 583   m_dot_air = v_dot_air * rho_air_manifold;
 584 }
 585
 586 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 587 /**
 588  * Calculate the fuel flow into the engine.
 589  *
 590  * Inputs: Mixture, thi_sea_level, p_amb_sea_level, p_amb, m_dot_air
 591  *
 592  * Outputs: equivalence_ratio, m_dot_fuel
 593  */
 594
 595 void FGPiston::doFuelFlow(void)
 596 {
 597   double thi_sea_level = 1.3 * Mixture; // Allows an AFR of infinity:1 to 11.3075:1
 598   equivalence_ratio = thi_sea_level; // * p_amb_sea_level / p_amb;
 599   double AFR = 10+(12*(1-Mixture));// mixture 10:1 to 22:1
 600   m_dot_fuel = m_dot_air / AFR;
 601   FuelFlow_gph = m_dot_fuel
 602     * 3600            // seconds to hours
 603     * 2.2046            // kg to lb
 604     / 6.0;            // lb to gal_us of gasoline
 605 //    / 6.6;            // lb to gal_us of kerosene
 606 }
 607
 608 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 609 /**
 610  * Calculate the power produced by the engine.
 611  *
 612  * Currently, the JSBSim propellor model does not allow the
 613  * engine to produce enough RPMs to get up to a high horsepower.
 614  * When tested with sufficient RPM, it has no trouble reaching
 615  * 200HP.
 616  *
 617  * Inputs: ManifoldPressure_inHg, p_amb, RPM, T_amb,
 618  *   Mixture_Efficiency_Correlation, Cycles, MaxHP
 619  *
 620  * Outputs: Percentage_Power, HP
 621  */
 622
 623 void FGPiston::doEnginePower(void)
 624 {
 625   if (Running) {
 626     double T_amb_degF = KelvinToFahrenheit(T_amb);
 627     double T_amb_sea_lev_degF = KelvinToFahrenheit(288);
 628
 629     // FIXME: this needs to be generalized
 630     double ME, Adjusted_BSFC;  // Convienience term for use in the calculations
 631     ME = Mixture_Efficiency_Correlation->GetValue(m_dot_fuel/m_dot_air);
 632     Adjusted_BSFC = (1/ThrottlePos) * BSFC;
 633     Percentage_Power = 1.000;
 634
 635     if ( Magnetos != 3 ) Percentage_Power *= SparkFailDrop;
 636
 637     HP = (FuelFlow_gph * 6.0 / Adjusted_BSFC )* ME * suction_loss * Percentage_Power;
 638
 639   } else {
 640
 641     // Power output when the engine is not running
 642     if (Cranking) {
 643       if (RPM < 10) {
 644         HP = StarterHP;
 645       } else if (RPM < IdleRPM*0.8) {
 646         HP = StarterHP + ((IdleRPM*0.8 - RPM) / 8.0);
 647         // This is a guess - would be nice to find a proper starter moter torque curve
 648       } else {
 649         HP = StarterHP;
 650       }
 651     } else {
 652       // Quick hack until we port the FMEP stuff
 653       if (RPM > 0.0)
 654         HP = -1.5;
 655       else
 656         HP = 0.0;
 657     }
 658   }
 659 //  cout << "Power = " << HP << "  RPM = " << RPM << "  Running = " << Running << "  Cranking = " << Cranking << endl;
 660 }
 661
 662 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 663 /**
 664  * Calculate the exhaust gas temperature.
 665  *
 666  * Inputs: equivalence_ratio, m_dot_fuel, calorific_value_fuel,
 667  *   Cp_air, m_dot_air, Cp_fuel, m_dot_fuel, T_amb, Percentage_Power
 668  *
 669  * Outputs: combustion_efficiency, ExhaustGasTemp_degK
 670  */
 671
 672 void FGPiston::doEGT(void)
 673 {
 674   double delta_T_exhaust;
 675   double enthalpy_exhaust;
 676   double heat_capacity_exhaust;
 677   double dEGTdt;
 678
 679   if ((Running) && (m_dot_air > 0.0)) {  // do the energy balance
 680     combustion_efficiency = Lookup_Combustion_Efficiency->GetValue(equivalence_ratio);
 681     enthalpy_exhaust = m_dot_fuel * calorific_value_fuel *
 682                               combustion_efficiency * 0.33;
 683     heat_capacity_exhaust = (Cp_air * m_dot_air) + (Cp_fuel * m_dot_fuel);
 684     delta_T_exhaust = enthalpy_exhaust / heat_capacity_exhaust;
 685     ExhaustGasTemp_degK = T_amb + delta_T_exhaust;
 686     ExhaustGasTemp_degK *= 0.444 + ((0.544 - 0.444) * Percentage_Power);
 687   } else {  // Drop towards ambient - guess an appropriate time constant for now
 688     combustion_efficiency = 0;
 689     dEGTdt = (RankineToKelvin(Atmosphere->GetTemperature()) - ExhaustGasTemp_degK) / 100.0;
 690     delta_T_exhaust = dEGTdt * dt;
 691     ExhaustGasTemp_degK += delta_T_exhaust;
 692   }
 693 }
 694
 695 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 696 /**
 697  * Calculate the cylinder head temperature.
 698  *
 699  * Inputs: T_amb, IAS, rho_air, m_dot_fuel, calorific_value_fuel,
 700  *   combustion_efficiency, RPM
 701  *
 702  * Outputs: CylinderHeadTemp_degK
 703  */
 704
 705 void FGPiston::doCHT(void)
 706 {
 707   double h1 = -95.0;
 708   double h2 = -3.95;
 709   double h3 = -0.05;
 710
 711   double arbitary_area = 1.0;
 712   double CpCylinderHead = 800.0;
 713   double MassCylinderHead = 8.0;
 714
 715   double temperature_difference = CylinderHeadTemp_degK - T_amb;
 716   double v_apparent = IAS * 0.5144444;
 717   double v_dot_cooling_air = arbitary_area * v_apparent;
 718   double m_dot_cooling_air = v_dot_cooling_air * rho_air;
 719   double dqdt_from_combustion =
 720     m_dot_fuel * calorific_value_fuel * combustion_efficiency * 0.33;
 721   double dqdt_forced = (h2 * m_dot_cooling_air * temperature_difference) +
 722     (h3 * RPM * temperature_difference);
 723   double dqdt_free = h1 * temperature_difference;
 724   double dqdt_cylinder_head = dqdt_from_combustion + dqdt_forced + dqdt_free;
 725
 726   double HeatCapacityCylinderHead = CpCylinderHead * MassCylinderHead;
 727
 728   CylinderHeadTemp_degK +=
 729     (dqdt_cylinder_head / HeatCapacityCylinderHead) * dt;
 730 }
 731
 732 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 733 /**
 734  * Calculate the oil temperature.
 735  *
 736  * Inputs: Percentage_Power, running flag.
 737  *
 738  * Outputs: OilTemp_degK
 739  */
 740
 741 void FGPiston::doOilTemperature(void)
 742 {
 743   double idle_percentage_power = 0.023;        // approximately
 744   double target_oil_temp;        // Steady state oil temp at the current engine conditions
 745   double time_constant;          // The time constant for the differential equation
 746
 747   if (Running) {
 748     target_oil_temp = 363;
 749     time_constant = 500;        // Time constant for engine-on idling.
 750     if (Percentage_Power > idle_percentage_power) {
 751       time_constant /= ((Percentage_Power / idle_percentage_power) / 10.0); // adjust for power
 752     }
 753   } else {
 754     target_oil_temp = RankineToKelvin(Atmosphere->GetTemperature());
 755     time_constant = 1000;  // Time constant for engine-off; reflects the fact
 756                            // that oil is no longer getting circulated
 757   }
 758
 759   double dOilTempdt = (target_oil_temp - OilTemp_degK) / time_constant;
 760
 761   OilTemp_degK += (dOilTempdt * dt);
 762 }
 763
 764 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 765 /**
 766  * Calculate the oil pressure.
 767  *
 768  * Inputs: RPM
 769  *
 770  * Outputs: OilPressure_psi
 771  */
 772
 773 void FGPiston::doOilPressure(void)
 774 {
 775   double Oil_Press_Relief_Valve = 60; // FIXME: may vary by engine
 776   double Oil_Press_RPM_Max = MaxRPM * 0.75;    // 75% of max rpm FIXME: may vary by engine
 777   double Design_Oil_Temp = 358;          // degK; FIXME: may vary by engine
 778   double Oil_Viscosity_Index = 0.25;
 779
 780   OilPressure_psi = (Oil_Press_Relief_Valve / Oil_Press_RPM_Max) * RPM;
 781
 782   if (OilPressure_psi >= Oil_Press_Relief_Valve) {
 783     OilPressure_psi = Oil_Press_Relief_Valve;
 784   }
 785
 786   OilPressure_psi += (Design_Oil_Temp - OilTemp_degK) * Oil_Viscosity_Index;
 787 }
 788
 789 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 790
 791 string FGPiston::GetEngineLabels(string delimeter)
 792 {
 793   std::ostringstream buf;
 794
 795   buf << Name << " Power Available (engine " << EngineNumber << " in HP)" << delimeter
 796       << Name << " HP (engine " << EngineNumber << ")" << delimeter
 797       << Name << " equivalent ratio (engine " << EngineNumber << ")" << delimeter
 798       << Name << " MAP (engine " << EngineNumber << ")" << delimeter
 799       << Thruster->GetThrusterLabels(EngineNumber, delimeter);
 800
 801   return buf.str();
 802 }
 803
 804 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 805
 806 string FGPiston::GetEngineValues(string delimeter)
 807 {
 808   std::ostringstream buf;
 809
 810   buf << PowerAvailable << delimeter << HP << delimeter
 811       << equivalence_ratio << delimeter << MAP << delimeter
 812       << Thruster->GetThrusterValues(EngineNumber, delimeter);
 813
 814   return buf.str();
 815 }
 816
 817 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 818 //
 819 //    The bitmasked value choices are as follows:
 820 //    unset: In this case (the default) JSBSim would only print
 821 //       out the normally expected messages, essentially echoing
 822 //       the config files as they are read. If the environment
 823 //       variable is not set, debug_lvl is set to 1 internally
 824 //    0: This requests JSBSim not to output any messages
 825 //       whatsoever.
 826 //    1: This value explicity requests the normal JSBSim
 827 //       startup messages
 828 //    2: This value asks for a message to be printed out when
 829 //       a class is instantiated
 830 //    4: When this value is set, a message is displayed when a
 831 //       FGModel object executes its Run() method
 832 //    8: When this value is set, various runtime state variables
 833 //       are printed out periodically
 834 //    16: When set various parameters are sanity checked and
 835 //       a message is printed out when they go out of bounds
 836
 837 void FGPiston::Debug(int from)
 838 {
 839   if (debug_lvl <= 0) return;
 840
 841   if (debug_lvl & 1) { // Standard console startup message output
 842     if (from == 0) { // Constructor
 843
 844       cout << "\n    Engine Name: "         << Name << endl;
 845       cout << "      MinManifoldPressure: " << MinManifoldPressure_inHg << endl;
 846       cout << "      MaxManifoldPressure: " << MaxManifoldPressure_inHg << endl;
 847       cout << "      MinMaP (Pa):         " << minMAP << endl;
 848       cout << "      MaxMaP (Pa): "         << maxMAP << endl;
 849       cout << "      Displacement: "        << Displacement             << endl;
 850       cout << "      MaxHP: "               << MaxHP                    << endl;
 851       cout << "      Cycles: "              << Cycles                   << endl;
 852       cout << "      IdleRPM: "             << IdleRPM                  << endl;
 853       cout << "      MaxThrottle: "         << MaxThrottle              << endl;
 854       cout << "      MinThrottle: "         << MinThrottle              << endl;
 855
 856       cout << endl;
 857       cout << "      Combustion Efficiency table:" << endl;
 858       Lookup_Combustion_Efficiency->Print();
 859       cout << endl;
 860
 861       cout << endl;
 862       cout << "      Power Mixture Correlation table:" << endl;
 863       Power_Mixture_Correlation->Print();
 864       cout << endl;
 865
 866       cout << endl;
 867       cout << "      Mixture Efficiency Correlation table:" << endl;
 868       Mixture_Efficiency_Correlation->Print();
 869       cout << endl;
 870
 871     }
 872   }
 873   if (debug_lvl & 2 ) { // Instantiation/Destruction notification
 874     if (from == 0) cout << "Instantiated: FGPiston" << endl;
 875     if (from == 1) cout << "Destroyed:    FGPiston" << endl;
 876   }
 877   if (debug_lvl & 4 ) { // Run() method entry print for FGModel-derived objects
 878   }
 879   if (debug_lvl & 8 ) { // Runtime state variables
 880   }
 881   if (debug_lvl & 16) { // Sanity checking
 882   }
 883   if (debug_lvl & 64) {
 884     if (from == 0) { // Constructor
 885       cout << IdSrc << endl;
 886       cout << IdHdr << endl;
 887     }
 888   }
 889 }
 890 } // namespace JSBSim