From 07f5955f8c01327a187eb86de7335b6b147d222a Mon Sep 17 00:00:00 2001 From: curt Date: Fri, 5 Oct 2001 20:27:01 +0000 Subject: [PATCH] Many engine model tweaks, updates, and enhancements. --- src/FDM/IO360.cxx | 937 +++++++++++++++++++++------------------------- src/FDM/IO360.hxx | 225 +++++------ 2 files changed, 523 insertions(+), 639 deletions(-) diff --git a/src/FDM/IO360.cxx b/src/FDM/IO360.cxx index f17b88e8e..b66e64386 100644 --- a/src/FDM/IO360.cxx +++ b/src/FDM/IO360.cxx @@ -1,10 +1,9 @@ -// Module: 10520c.c -// Author: Phil Schubert -// Date started: 12/03/99 -// Purpose: Models a Continental IO-520-M Engine -// Called by: FGSimExec -// -// Copyright (C) 1999 Philip L. Schubert (philings@ozemail.com.au) +// IO360.cxx - a piston engine model currently for the IO360 engine fitted to the C172 +// but with the potential to model other naturally aspirated piston engines +// given appropriate config input. +// +// Written by David Luff, started 2000. +// Based on code by Phil Schubert, started 1999. // // This program is free software; you can redistribute it and/or // modify it under the terms of the GNU General Public License as @@ -18,79 +17,7 @@ // // You should have received a copy of the GNU General Public License // along with this program; if not, write to the Free Software -// Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA -// 02111-1307, USA. -// -// Further information about the GNU General Public License can also -// be found on the world wide web at http://www.gnu.org. -// -// FUNCTIONAL DESCRIPTION -// ------------------------------------------------------------------------ -// Models a Continental IO-520-M engine. This engine is used in Cessna -// 210, 310, Beechcraft Bonaza and Baron C55. The equations used below -// were determined by a first and second order curve fits using Excel. -// The data is from the Cessna Aircraft Corporations Engine and Flight -// Computer for C310. Part Number D3500-13 -// -// ARGUMENTS -// ------------------------------------------------------------------------ -// -// -// HISTORY -// ------------------------------------------------------------------------ -// 12/03/99 PLS Created -// 07/03/99 PLS Added Calculation of Density, and Prop_Torque -// 07/03/99 PLS Restructered Variables to allow easier implementation -// of Classes -// 15/03/99 PLS Added Oil Pressure, Oil Temperature and CH Temp -// ------------------------------------------------------------------------ -// INCLUDES -// ------------------------------------------------------------------------ -// -// -///////////////////////////////////////////////////////////////////// -// -// Modified by Dave Luff (david.luff@nottingham.ac.uk) September 2000 -// -// Altered manifold pressure range to add a minimum value at idle to simulate the throttle stop / idle bypass valve, -// and to reduce the maximum value whilst the engine is running to slightly below ambient to account for CdA losses across the throttle -// -// Altered it a bit to model an IO360 from C172 - 360 cubic inches, 180 HP max, fixed pitch prop -// Added a simple fixed pitch prop model by Nev Harbor - this is not intended as a final model but simply a hack to get it running for now -// I used Phil's ManXRPM correlation for power rather than do a new one for the C172 for now, but altered it a bit to reduce power at the low end -// -// Added EGT model based on combustion efficiency and an energy balance with the exhaust gases -// -// Added a mixture - power correlation based on a curve in the IO360 operating manual -// -// I've tried to match the prop and engine model to give roughly 600 RPM idle and 180 HP at 2700 RPM -// but it is by no means currently at a completed stage - DCL 15/9/00 -// -// DCL 28/09/00 - Added estimate of engine and prop inertia and changed engine speed calculation to be calculated from Angular acceleration = Torque / Inertia. -// Requires a timestep to be passed to FGNewEngine::init and currently assumes this timestep does not change. -// Could easily be altered to pass a variable timestep to FGNewEngine::update every step instead if required. -// -// DCL 27/10/00 - Added first stab at cylinder head temperature model -// See the comment block in the code for details -// -// DCL 02/11/00 - Modified EGT code to reduce values to those more representative of a sensor downstream -// -// DCL 02/02/01 - Changed the prop model to one based on efficiency and co-efficient of power curves from McCormick instead of the -// blade element method we were using previously. This works much better, and is similar to how Jon is doing it in JSBSim. -// -// DCL 08/02/01 - Overhauled fuel consumption rate support. -// -// DCL 22/03/01 - Added input of actual air pressure and temperature (and hence density) to the model. Hence the power correlation -// with pressure height and temperature is no longer required since the power is based on the actual manifold pressure. -// -// DCL 22/03/01 - based on Riley's post on the list (25 rpm gain at 1000 rpm as lever is pulled out from full rich) -// I have reduced the sea level full rich mixture to thi = 1.3 -// -// DCL 18/9/01 - Got the engine to start and stop in response to the magneto switch. -// Changed all PI to LS_PI (in ls_constants.h). -// Engine now checks for fuel and stops when not available. -// -////////////////////////////////////////////////////////////////////// +// Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. #include @@ -108,46 +35,288 @@ SG_USING_STD(cout); #include
-// Static utility functions +//************************************************************************************* +// Initialise the engine model +void FGNewEngine::init(double dt) { -// Calculate Density Ratio -static float Density_Ratio ( float x ) -{ - float y ; - y = ((3E-10 * x * x) - (3E-05 * x) + 0.9998); - return(y); + // These constants should probably be moved eventually + CONVERT_CUBIC_INCHES_TO_METERS_CUBED = 1.638706e-5; + CONVERT_HP_TO_WATTS = 745.6999; + + // Properties of working fluids + Cp_air = 1005; // J/KgK + Cp_fuel = 1700; // J/KgK + calorific_value_fuel = 47.3e6; // W/Kg Note that this is only an approximate value + rho_fuel = 800; // kg/m^3 - an estimate for now + R_air = 287.3; + + // environment inputs + p_amb_sea_level = 101325; // Pascals + + // Control inputs - ARE THESE NEEDED HERE??? + Throttle_Lever_Pos = 75; + Propeller_Lever_Pos = 75; + Mixture_Lever_Pos = 100; + + //misc + IAS = 0; + time_step = dt; + + // Engine Specific Variables that should be read in from a config file + MaxHP = 200; //Lycoming IO360 -A-C-D series +// MaxHP = 180; //Current Lycoming IO360 ? +// displacement = 520; //Continental IO520-M + displacement = 360; //Lycoming IO360 + displacement_SI = displacement * CONVERT_CUBIC_INCHES_TO_METERS_CUBED; + engine_inertia = 0.2; //kgm^2 - value taken from a popular family saloon car engine - need to find an aeroengine value !!!!! + prop_inertia = 0.05; //kgm^2 - this value is a total guess - dcl + Max_Fuel_Flow = 130; // Units??? Do we need this variable any more?? + + // Engine specific variables that maybe should be read in from config but are pretty generic and won't vary much for a naturally aspirated piston engine. + Max_Manifold_Pressure = 28.50; //Inches Hg. An approximation - should be able to find it in the engine performance data + Min_Manifold_Pressure = 6.5; //Inches Hg. This is a guess corresponding to approx 0.24 bar MAP (7 in Hg) - need to find some proper data for this + Max_RPM = 2700; + Min_RPM = 600; //Recommended idle from Continental data sheet + Mag_Derate_Percent = 5; + Gear_Ratio = 1; + n_R = 2; // Number of crank revolutions per power cycle - 2 for a 4 stroke engine. + + // Various bits of housekeeping describing the engines initial state. + running = fgGetBool("/engines/engine[0]/running"); + cranking = false; + crank_counter = false; + fgSetBool("/engines/engine[0]/cranking", false); + + // Initialise Engine Variables used by this instance + if(running) + RPM = 600; + else + RPM = 0; + Percentage_Power = 0; + Manifold_Pressure = 29.96; // Inches + Fuel_Flow_gals_hr = 0; +// Torque = 0; + Torque_SI = 0; + CHT_degK = 298.0; //deg Kelvin + CHT_degF = (CHT_degF * 1.8) - 459.67; //deg Fahrenheit + Mixture = 14; + Oil_Pressure = 0; // PSI + Oil_Temp = 85; // Deg C + current_oil_temp = 298.0; //deg Kelvin + /**** one of these is superfluous !!!!***/ + HP = 0; + RPS = 0; + Torque_Imbalance = 0; + + // Initialise Propellor Variables used by this instance + FGProp1_RPS = 0; + // Hardcode propellor for now - the following two should be read from config eventually + prop_diameter = 1.8; // meters + blade_angle = 23.0; // degrees } +//***************************************************************************** +// update the engine model based on current control positions +void FGNewEngine::update() { -// Calculate Air Density - Rho, using the ideal gas equation -// Takes and returns SI values -static float Density ( float temperature, float pressure ) -{ - // rho = P / RT - // R = 287.3 for air +/* + // Hack for testing - should output every 5 seconds + static int count1 = 0; + if(count1 == 0) { +// cout << "P_atmos = " << p_amb << " T_atmos = " << T_amb << '\n'; +// cout << "Manifold pressure = " << Manifold_Pressure << " True_Manifold_Pressure = " << True_Manifold_Pressure << '\n'; +// cout << "p_amb_sea_level = " << p_amb_sea_level << '\n'; +// cout << "equivalence_ratio = " << equivalence_ratio << '\n'; +// cout << "combustion_efficiency = " << combustion_efficiency << '\n'; +// cout << "AFR = " << 14.7 / equivalence_ratio << '\n'; +// cout << "Mixture lever = " << Mixture_Lever_Pos << '\n'; +// cout << "n = " << RPM << " rpm\n"; +// cout << "T_amb = " << T_amb << '\n'; +// cout << "running = " << running << '\n'; +// cout << "fuel = " << fgGetFloat("/consumables/fuel/tank[0]/level-gal_us") << '\n'; + cout << "Percentage_Power = " << Percentage_Power << '\n'; + cout << "current_oil_temp = " << current_oil_temp << '\n'; + } + count1++; + if(count1 == 600) + count1 = 0; +*/ - float R = 287.3; - float rho = pressure / (R * temperature); - return(rho); + // Check parameters that may alter the operating state of the engine. + // (spark, fuel, starter motor etc) + + // Check for spark + bool Magneto_Left = false; + bool Magneto_Right = false; + int mag_pos = fgGetInt("/engines/engine[0]/magneto"); + // Magneto positions: + // 0 -> off + // 1 -> left only + // 2 -> right only + // 3 -> both + if(mag_pos != 0) { + spark = true; + } else { + spark = false; + } // neglects battery voltage, master on switch, etc for now. + if((mag_pos == 1) || (mag_pos > 2)) + Magneto_Left = true; + if(mag_pos > 1) + Magneto_Right = true; + + // crude check for fuel + if((fgGetFloat("/consumables/fuel/tank[0]/level-gal_us") > 0) || (fgGetFloat("/consumables/fuel/tank[1]/level-gal_us") > 0)) { + fuel = true; + } else { + fuel = false; + } // Need to make this better, eg position of fuel selector switch. + + // Check if we are turning the starter motor + bool temp = fgGetBool("/engines/engine[0]/starter"); + if(cranking != temp) { + // This check saves .../cranking from getting updated every loop - they only update when changed. + cranking = temp; + if(cranking) + fgSetBool("/engines/engine[0]/cranking", true); + else + fgSetBool("/engines/engine[0]/cranking", false); + crank_counter = 0; + } + // Note that although /engines/engine[0]/starter and /engines/engine[0]/cranking might appear to be duplication it is + // not since the starter may be engaged with the battery voltage too low for cranking to occur (or perhaps the master + // switch just left off) and the sound manager will read .../cranking to determine wether to play a cranking sound. + // For now though none of that is implemented so cranking can be set equal to .../starter without further checks. + +// int Alternate_Air_Pos =0; // Off = 0. Reduces power by 3 % for same throttle setting + // DCL - don't know what this Alternate_Air_Pos is - this is a leftover from the Schubert code. + + //Check mode of engine operation + if(cranking) { + crank_counter++; + if(RPM <= 480) { + RPM += 100; + if(RPM > 480) + RPM = 480; + } else { + // consider making a horrible noise if the starter is engaged with the engine running + } + } + if((!running) && (spark) && (fuel) && (crank_counter > 120)) { + // start the engine if revs high enough + if(RPM > 450) { + // For now just instantaneously start but later we should maybe crank for a bit + running = true; + fgSetBool("/engines/engine[0]/running", true); +// RPM = 600; + } + } + if( (running) && ((!spark)||(!fuel)) ) { + // Cut the engine + // note that we only cut the power - the engine may continue to spin if the prop is in a moving airstream + running = false; + fgSetBool("/engines/engine[0]/running", false); + } + + // Now we've ascertained whether the engine is running or not we can start to do the engine calculations 'proper' + + // Calculate Sea Level Manifold Pressure + Manifold_Pressure = Calc_Manifold_Pressure( Throttle_Lever_Pos, Max_Manifold_Pressure, Min_Manifold_Pressure ); + // cout << "manifold pressure = " << Manifold_Pressure << endl; + + //Then find the actual manifold pressure (the calculated one is the sea level pressure) + True_Manifold_Pressure = Manifold_Pressure * p_amb / p_amb_sea_level; + + //Do the fuel flow calculations + Calc_Fuel_Flow_Gals_Hr(); + + //Calculate engine power + Calc_Percentage_Power(Magneto_Left, Magneto_Right); + HP = Percentage_Power * MaxHP / 100.0; + Power_SI = HP * CONVERT_HP_TO_WATTS; + + // FMEP calculation. For now we will just use this during motored operation. + // Eventually we will calculate IMEP and use the FMEP all the time to give BMEP (maybe!) + if(!running) { + // This FMEP data is from the Patton paper, assumes fully warm conditions. + FMEP = 1e-12*pow(RPM,4) - 1e-8*pow(RPM,3) + 5e-5*pow(RPM,2) - 0.0722*RPM + 154.85; + // Gives FMEP in kPa - now convert to Pa + FMEP *= 1000; + } else { + FMEP = 0.0; + } + // Is this total FMEP or friction FMEP ??? + + Torque_FMEP = (FMEP * displacement_SI) / (2.0 * LS_PI * n_R); + + // Calculate Engine Torque. Check for div by zero since percentage power correlation does not guarantee zero power at zero rpm. + // However this is problematical since there is a resistance to movement even at rest + // Ie this is a dynamics equation not a statics one. This can be solved by going over to MEP based torque calculations. + if(RPM == 0) { + Torque_SI = 0 - Torque_FMEP; + } + else { + Torque_SI = ((Power_SI * 60.0) / (2.0 * LS_PI * RPM)) - Torque_FMEP; //Torque = power / angular velocity + // cout << Torque << " Nm\n"; + } + + //Calculate Exhaust gas temperature + Calc_EGT(); + + // Calculate Cylinder Head Temperature + CHT_degK = Calc_CHT(CHT_degK); + CHT_degF = (CHT_degK * 1.8) - 459.67; + + // Calculate oil temperature + current_oil_temp = Calc_Oil_Temp(current_oil_temp); + + // Calculate Oil Pressure + Oil_Pressure = Calc_Oil_Press( Oil_Temp, RPM ); + + // Now do the Propeller Calculations + Do_Prop_Calcs(); + +// Now do the engine - prop torque balance to calculate final RPM + + //Calculate new RPM from torque balance and inertia. + Torque_Imbalance = Torque_SI - prop_torque; //This gives a +ve value when the engine torque exeeds the prop torque + // (Engine torque is +ve when it acts in the direction of engine revolution, prop torque is +ve when it opposes the direction of engine revolution) + + angular_acceleration = Torque_Imbalance / (engine_inertia + prop_inertia); + angular_velocity_SI += (angular_acceleration * time_step); + // Don't let the engine go into reverse + if(angular_velocity_SI < 0) + angular_velocity_SI = 0; + RPM = (angular_velocity_SI * 60) / (2.0 * LS_PI); + + // And finally a last check on the engine state after doing the torque balance with the prop - have we stalled? + if(running) { + //Check if we have stalled the engine + if (RPM == 0) { + running = false; + fgSetBool("/engines/engine[0]/running", false); + } else if((RPM <= 480) && (cranking)) { + //Make sure the engine noise dosn't play if the engine won't start due to eg mixture lever pulled out. + running = false; + fgSetBool("/engines/engine[0]/running", false); + } + } } - -// Calculate Speed in FPS given Knots CAS -static float IAS_to_FPS (float x) -{ - float y; - y = x * 1.68888888; - return y; -} +//***************************************************************************************************** // FGNewEngine member functions +//////////////////////////////////////////////////////////////////////////////////////////// +// Return the combustion efficiency as a function of equivalence ratio +// +// Combustion efficiency values from Heywood, +// "Internal Combustion Engine Fundamentals", ISBN 0-07-100499-8 +//////////////////////////////////////////////////////////////////////////////////////////// float FGNewEngine::Lookup_Combustion_Efficiency(float thi_actual) { const int NUM_ELEMENTS = 11; float thi[NUM_ELEMENTS] = {0.0, 0.9, 1.0, 1.05, 1.1, 1.15, 1.2, 1.3, 1.4, 1.5, 1.6}; //array of equivalence ratio values float neta_comb[NUM_ELEMENTS] = {0.98, 0.98, 0.97, 0.95, 0.9, 0.85, 0.79, 0.7, 0.63, 0.57, 0.525}; //corresponding array of combustion efficiency values - //combustion efficiency values from Heywood, "Internal Combustion Engine Fundamentals", ISBN 0-07-100499-8 float neta_comb_actual = 0.0f; float factor; @@ -233,15 +402,86 @@ float FGNewEngine::Power_Mixture_Correlation(float thi_actual) return mixPerPow_actual; } +// Calculate Cylinder Head Temperature +// Crudely models the cylinder head as an arbitary lump of arbitary size and area with one third of combustion energy +// as heat input and heat output as a function of airspeed and temperature. Could be improved!!! +float FGNewEngine::Calc_CHT(float CHT) +{ + float h1 = -95.0; //co-efficient for free convection + float h2 = -3.95; //co-efficient for forced convection + float h3 = -0.05; //co-efficient for forced convection due to prop backwash + float v_apparent; //air velocity over cylinder head in m/s + float v_dot_cooling_air; + float m_dot_cooling_air; + float temperature_difference; + float arbitary_area = 1.0; + float dqdt_from_combustion; + float dqdt_forced; //Rate of energy transfer to/from cylinder head due to forced convection (Joules) (sign convention: to cylinder head is +ve) + float dqdt_free; //Rate of energy transfer to/from cylinder head due to free convection (Joules) (sign convention: to cylinder head is +ve) + float dqdt_cylinder_head; //Overall energy change in cylinder head + float CpCylinderHead = 800.0; //FIXME - this is a guess - I need to look up the correct value + float MassCylinderHead = 8.0; //Kg - this is a guess - it dosn't have to be absolutely accurate but can be varied to alter the warm-up rate + float HeatCapacityCylinderHead; + float dCHTdt; + // The above values are hardwired to give reasonable results for an IO360 (C172 engine) + // Now adjust to try to compensate for arbitary sized engines - this is currently untested + arbitary_area *= (MaxHP / 180.0); + MassCylinderHead *= (MaxHP / 180.0); + + temperature_difference = CHT - T_amb; + + v_apparent = IAS * 0.5144444; //convert from knots to m/s + v_dot_cooling_air = arbitary_area * v_apparent; + m_dot_cooling_air = v_dot_cooling_air * rho_air; + + //Calculate rate of energy transfer to cylinder head from combustion + dqdt_from_combustion = m_dot_fuel * calorific_value_fuel * combustion_efficiency * 0.33; + + //Calculate rate of energy transfer from cylinder head due to cooling NOTE is calculated as rate to but negative + dqdt_forced = (h2 * m_dot_cooling_air * temperature_difference) + (h3 * RPM * temperature_difference); + dqdt_free = h1 * temperature_difference; + + //Calculate net rate of energy transfer to or from cylinder head + dqdt_cylinder_head = dqdt_from_combustion + dqdt_forced + dqdt_free; + + HeatCapacityCylinderHead = CpCylinderHead * MassCylinderHead; + + dCHTdt = dqdt_cylinder_head / HeatCapacityCylinderHead; + + CHT += (dCHTdt * time_step); + + return(CHT); +} + +// Calculate exhaust gas temperature +void FGNewEngine::Calc_EGT() +{ + combustion_efficiency = Lookup_Combustion_Efficiency(equivalence_ratio); //The combustion efficiency basically tells us what proportion of the fuels calorific value is released + + //now calculate energy release to exhaust + //We will assume a three way split of fuel energy between useful work, the coolant system and the exhaust system + //This is a reasonable first suck of the thumb estimate for a water cooled automotive engine - whether it holds for an air cooled aero engine is probably open to question + //Regardless - it won't affect the variation of EGT with mixture, and we can always put a multiplier on EGT to get a reasonable peak value. + enthalpy_exhaust = m_dot_fuel * calorific_value_fuel * combustion_efficiency * 0.33; + heat_capacity_exhaust = (Cp_air * m_dot_air) + (Cp_fuel * m_dot_fuel); + delta_T_exhaust = enthalpy_exhaust / heat_capacity_exhaust; + + EGT = T_amb + delta_T_exhaust; + + //The above gives the exhaust temperature immediately prior to leaving the combustion chamber + //Now derate to give a more realistic figure as measured downstream + //For now we will aim for a peak of around 400 degC (750 degF) + + EGT *= 0.444 + ((0.544 - 0.444) * Percentage_Power / 100.0); + + EGT_degF = (EGT * 1.8) - 459.67; +} // Calculate Manifold Pressure based on Throttle lever Position float FGNewEngine::Calc_Manifold_Pressure ( float LeverPosn, float MaxMan, float MinMan) { float Inches; - // if ( x < = 0 ) { - // x = 0.00001; - // } //Note that setting the manifold pressure as a function of lever position only is not strictly accurate //MAP is also a function of engine speed. (and ambient pressure if we are going for an actual MAP model) @@ -251,232 +491,16 @@ float FGNewEngine::Calc_Manifold_Pressure ( float LeverPosn, float MaxMan, float if(Inches < MinMan) Inches = MinMan; + //Check for stopped engine (crudest way of compensating for engine speed) + if(RPM == 0) + Inches = 29.92; + return Inches; } - - - -// Calculate Oil Temperature -float FGNewEngine::Calc_Oil_Temp (float Fuel_Flow, float Mixture, float IAS) +// Calculate fuel flow in gals/hr +void FGNewEngine::Calc_Fuel_Flow_Gals_Hr() { - float Oil_Temp = 85; - - return (Oil_Temp); -} - -// Calculate Oil Pressure -float FGNewEngine::Calc_Oil_Press (float Oil_Temp, float Engine_RPM) -{ - float Oil_Pressure = 0; //PSI - float Oil_Press_Relief_Valve = 60; //PSI - float Oil_Press_RPM_Max = 1800; - float Design_Oil_Temp = 85; //Celsius - float Oil_Viscosity_Index = 0.25; // PSI/Deg C - float Temp_Deviation = 0; // Deg C - - Oil_Pressure = (Oil_Press_Relief_Valve / Oil_Press_RPM_Max) * Engine_RPM; - - // Pressure relief valve opens at Oil_Press_Relief_Valve PSI setting - if (Oil_Pressure >= Oil_Press_Relief_Valve) - { - Oil_Pressure = Oil_Press_Relief_Valve; - } - - // Now adjust pressure according to Temp which affects the viscosity - - Oil_Pressure += (Design_Oil_Temp - Oil_Temp) * Oil_Viscosity_Index; - - return Oil_Pressure; -} - -//************************************************************************************* -// Initialise the engine model -void FGNewEngine::init(double dt) { - - // These constants should probably be moved eventually - CONVERT_CUBIC_INCHES_TO_METERS_CUBED = 1.638706e-5; - CONVERT_HP_TO_WATTS = 745.6999; - - // Properties of working fluids - Cp_air = 1005; // J/KgK - Cp_fuel = 1700; // J/KgK - calorific_value_fuel = 47.3e6; // W/Kg Note that this is only an approximate value - rho_fuel = 800; // kg/m^3 - an estimate for now - R_air = 287.3; - p_amb_sea_level = 101325; - - // Control and environment inputs - IAS = 0; - Throttle_Lever_Pos = 75; - Propeller_Lever_Pos = 75; - Mixture_Lever_Pos = 100; - - time_step = dt; - - // Engine Specific Variables. - // Will be set in a parameter file to be read in to create - // and instance for each engine. - Max_Manifold_Pressure = 28.50; //Inches Hg. An approximation - should be able to find it in the engine performance data - Min_Manifold_Pressure = 6.5; //Inches Hg. This is a guess corresponding to approx 0.24 bar MAP (7 in Hg) - need to find some proper data for this - Max_RPM = 2700; - Min_RPM = 600; //Recommended idle from Continental data sheet - Max_Fuel_Flow = 130; - Mag_Derate_Percent = 5; -// MaxHP = 285; //Continental IO520-M - MaxHP = 200; //Lycoming IO360 -A-C-D series -// MaxHP = 180; //Current Lycoming IO360 ? -// displacement = 520; //Continental IO520-M - displacement = 360; //Lycoming IO360 - displacement_SI = displacement * CONVERT_CUBIC_INCHES_TO_METERS_CUBED; - engine_inertia = 0.2; //kgm^2 - value taken from a popular family saloon car engine - need to find an aeroengine value !!!!! - prop_inertia = 0.05; //kgm^2 - this value is a total guess - dcl - Gear_Ratio = 1; - n_R = 2; // Number of crank revolutions per power cycle - 2 for a 4 cylinder engine. - - running = fgGetBool("/engines/engine[0]/running"); - cranking = false; - fgSetBool("/engines/engine[0]/cranking", false); - - // Initialise Engine Variables used by this instance - if(running) - RPM = 600; - else - RPM = 0; - Percentage_Power = 0; - Manifold_Pressure = 29.00; // Inches - Fuel_Flow_gals_hr = 0; - Torque = 0; - Torque_SI = 0; - CHT = 298.0; //deg Kelvin - CHT_degF = (CHT * 1.8) - 459.67; //deg Fahrenheit - Mixture = 14; - Oil_Pressure = 0; // PSI - Oil_Temp = 85; // Deg C - HP = 0; - RPS = 0; - Torque_Imbalance = 0; - - // Initialise Propellor Variables used by this instance - FGProp1_Thrust = 0; - FGProp1_RPS = 0; - FGProp1_Blade_Angle = 13.5; - prop_diameter = 1.8; // meters - blade_angle = 23.0; // degrees -} - - -//***************************************************************************** -//***************************************************************************** -// update the engine model based on current control positions -void FGNewEngine::update() { - -/* - // Hack for testing - should output every 5 seconds - static int count1 = 0; - if(count1 == 0) { -// cout << "P_atmos = " << p_amb << " T_atmos = " << T_amb << '\n'; -// cout << "Manifold pressure = " << Manifold_Pressure << " True_Manifold_Pressure = " << True_Manifold_Pressure << '\n'; -// cout << "p_amb_sea_level = " << p_amb_sea_level << '\n'; -// cout << "equivalence_ratio = " << equivalence_ratio << '\n'; -// cout << "combustion_efficiency = " << combustion_efficiency << '\n'; -// cout << "AFR = " << 14.7 / equivalence_ratio << '\n'; -// cout << "Mixture lever = " << Mixture_Lever_Pos << '\n'; -// cout << "n = " << RPM << " rpm\n"; -// cout << "T_amb = " << T_amb << '\n'; - cout << "running = " << running << '\n'; - cout << "fuel = " << fgGetFloat("/consumables/fuel/tank[0]/level-gal_us") << '\n'; - } - count1++; - if(count1 == 600) - count1 = 0; -*/ - - float ManXRPM = 0; - float Vo = 0; - float V1 = 0; - - // Parameters that alter the operation of the engine. (spark, fuel, starter motor etc) - // Check for spark - int Magneto_Left = 0; - int Magneto_Right = 0; - int mag_pos = fgGetInt("/engines/engine[0]/magneto"); - // Magneto positions: - // 0 -> off - // 1 -> left only - // 2 -> right only - // 3 -> both - if(mag_pos != 0) { - spark = true; - } else { - spark = false; - } // neglects battery voltage, master on switch, etc for now. - if((mag_pos == 1) || (mag_pos > 2)) - Magneto_Left = 1; - if(mag_pos > 1) - Magneto_Right = 1; - - // crude check for fuel - if((fgGetFloat("/consumables/fuel/tank[0]/level-gal_us") > 0) || (fgGetFloat("/consumables/fuel/tank[1]/level-gal_us") > 0)) { - fuel = true; - } else { - fuel = false; - } // Need to make this better, eg position of fuel selector switch. - - // Check if we are turning the starter motor - bool temp = fgGetBool("/engines/engine[0]/starter"); - if(cranking != temp) { - // This check saves .../cranking from getting updated every loop - they only update when changed. - cranking = temp; - if(cranking) - fgSetBool("/engines/engine[0]/cranking", true); - else - fgSetBool("/engines/engine[0]/cranking", false); - } - // Note that although /engines/engine[0]/starter and /engines/engine[0]/cranking might appear to be duplication it is - // not since the starter may be engaged with the battery voltage too low for cranking to occur (or perhaps the master - // switch just left off) and the sound manager will read .../cranking to determine wether to play a cranking sound. - // For now though none of that is implemented so cranking can be set equal to .../starter without further checks. - - int Alternate_Air_Pos =0; // Off = 0. Reduces power by 3 % for same throttle setting - // DCL - don't know what this Alternate_Air_Pos is - this is a leftover from the Schubert code. - - //Check mode of engine operation - if(cranking) { - if(RPM <= 480) { - RPM += 100; - if(RPM > 480) - RPM = 480; - } else { - // consider making a horrible noise if the starter is engaged with the engine running - } - } - if((!running) && (spark) && (fuel)) { - // start the engine if revs high enough - if(RPM > 450) { - // For now just instantaneously start but later we should maybe crank for a bit - running = true; - fgSetBool("/engines/engine[0]/running", true); - RPM = 600; - } - } - if( (running) && ((!spark)||(!fuel)) ) { - // Cut the engine - // note that we only cut the power - the engine may continue to spin if the prop is in a moving airstream - running = false; - fgSetBool("/engines/engine[0]/running", false); - } - - // Calculate Sea Level Manifold Pressure - Manifold_Pressure = Calc_Manifold_Pressure( Throttle_Lever_Pos, Max_Manifold_Pressure, Min_Manifold_Pressure ); - // cout << "manifold pressure = " << Manifold_Pressure << endl; - - //Then find the actual manifold pressure (the calculated one is the sea level pressure) - True_Manifold_Pressure = Manifold_Pressure * p_amb / p_amb_sea_level; - -//************* -//DCL - next calculate m_dot_air and m_dot_fuel into engine - //DCL - calculate mass air flow into engine based on speed and load - separate this out into a function eventually //t_amb is actual temperature calculated from altitude //calculate density from ideal gas equation @@ -501,16 +525,16 @@ void FGNewEngine::update() { equivalence_ratio = thi_sea_level * p_amb_sea_level / p_amb; //ie as we go higher the mixture gets richer for a given lever position m_dot_fuel = m_dot_air / 14.7 * equivalence_ratio; Fuel_Flow_gals_hr = (m_dot_fuel / rho_fuel) * 264.172 * 3600.0; // Note this assumes US gallons +} -//*********************************************************************** -//Engine power and torque calculations - +// Calculate the percentage of maximum rated power delivered as a function of Manifold pressure multiplied by engine speed (rpm) +// This is not necessarilly the best approach but seems to work for now. +// May well need tweaking at the bottom end if the prop model is changed. +void FGNewEngine::Calc_Percentage_Power(bool mag_left, bool mag_right) +{ // For a given Manifold Pressure and RPM calculate the % Power // Multiply Manifold Pressure by RPM - ManXRPM = True_Manifold_Pressure * RPM; -// ManXRPM = Manifold_Pressure * RPM; - // cout << ManXRPM; - // cout << endl; + float ManXRPM = True_Manifold_Pressure * RPM; /* // Phil's %power correlation @@ -549,193 +573,92 @@ void FGNewEngine::update() { if(!running) { // cout << "Both OFF\n"; Percentage_Power = 0; - } else if (Magneto_Left && Magneto_Right) { + } else if (mag_left && mag_right) { // cout << "Both On "; - } else if (Magneto_Left == 0 || Magneto_Right== 0) { + } else if (mag_left == 0 || mag_right== 0) { // cout << "1 Magneto Failed "; Percentage_Power = Percentage_Power * ((100.0 - Mag_Derate_Percent)/100.0); // cout << FGEng1_Percentage_Power << "%" << "\t"; } - +/* //DCL - stall the engine if RPM drops below 450 - this is possible if mixture lever is pulled right out //This is a kludge that I should eliminate by adding a total fmep estimation. if(RPM < 450) Percentage_Power = 0; - +*/ if(Percentage_Power < 0) Percentage_Power = 0; +} - // FMEP calculation. For now we will just use this during motored operation, ie when %Power == 0. - // Eventually we will calculate IMEP and use the FMEP all the time to give BMEP - // - if(Percentage_Power == 0) { - // This FMEP data is from the Patton paper, assumes fully warm conditions. - FMEP = 1e-12*pow(RPM,4) - 1e-8*pow(RPM,3) + 5e-5*pow(RPM,2) - 0.0722*RPM + 154.85; - // Gives FMEP in kPa - now convert to Pa - FMEP *= 1000; +// Calculate Oil Temperature in degrees Kelvin +float FGNewEngine::Calc_Oil_Temp (float oil_temp) +{ + float idle_percentage_power = 2.3; // approximately + float target_oil_temp; // Steady state oil temp at the current engine conditions + float time_constant; // The time constant for the differential equation + if(running) { + target_oil_temp = 363; + time_constant = 500; // Time constant for engine-on idling. + if(Percentage_Power > idle_percentage_power) { + time_constant /= ((Percentage_Power / idle_percentage_power) / 10.0); // adjust for power + } } else { - FMEP = 0.0; + target_oil_temp = 298; + time_constant = 1000; // Time constant for engine-off; reflects the fact that oil is no longer getting circulated } - Torque_FMEP = (FMEP * displacement_SI) / (2.0 * LS_PI * n_R); + float dOilTempdt = (target_oil_temp - oil_temp) / time_constant; - HP = Percentage_Power * MaxHP / 100.0; + oil_temp += (dOilTempdt * time_step); - Power_SI = HP * CONVERT_HP_TO_WATTS; + return (oil_temp); +} - // Calculate Engine Torque. Check for div by zero since percentage power correlation does not guarantee zero power at zero rpm. - // However this is problematical since there is a resistance to movement even at rest - // Ie this is a dynamics equation not a statics one. This can be solved by going over to MEP based torque calculations. - if(RPM == 0) { - Torque_SI = 0 - Torque_FMEP; - } - else { - Torque_SI = ((Power_SI * 60.0) / (2.0 * LS_PI * RPM)) - Torque_FMEP; //Torque = power / angular velocity - // cout << Torque << " Nm\n"; +// Calculate Oil Pressure +float FGNewEngine::Calc_Oil_Press (float Oil_Temp, float Engine_RPM) +{ + float Oil_Pressure = 0; //PSI + float Oil_Press_Relief_Valve = 60; //PSI + float Oil_Press_RPM_Max = 1800; + float Design_Oil_Temp = 85; //Celsius + float Oil_Viscosity_Index = 0.25; // PSI/Deg C +// float Temp_Deviation = 0; // Deg C + + Oil_Pressure = (Oil_Press_Relief_Valve / Oil_Press_RPM_Max) * Engine_RPM; + + // Pressure relief valve opens at Oil_Press_Relief_Valve PSI setting + if (Oil_Pressure >= Oil_Press_Relief_Valve) { + Oil_Pressure = Oil_Press_Relief_Valve; } -//********************************************************************** -//Calculate Exhaust gas temperature + // Now adjust pressure according to Temp which affects the viscosity - // cout << "Thi = " << equivalence_ratio << '\n'; + Oil_Pressure += (Design_Oil_Temp - Oil_Temp) * Oil_Viscosity_Index; - combustion_efficiency = Lookup_Combustion_Efficiency(equivalence_ratio); //The combustion efficiency basically tells us what proportion of the fuels calorific value is released + return Oil_Pressure; +} - // cout << "Combustion efficiency = " << combustion_efficiency << '\n'; - //now calculate energy release to exhaust - //We will assume a three way split of fuel energy between useful work, the coolant system and the exhaust system - //This is a reasonable first suck of the thumb estimate for a water cooled automotive engine - whether it holds for an air cooled aero engine is probably open to question - //Regardless - it won't affect the variation of EGT with mixture, and we can always put a multiplier on EGT to get a reasonable peak value. - enthalpy_exhaust = m_dot_fuel * calorific_value_fuel * combustion_efficiency * 0.33; - heat_capacity_exhaust = (Cp_air * m_dot_air) + (Cp_fuel * m_dot_fuel); - delta_T_exhaust = enthalpy_exhaust / heat_capacity_exhaust; -// delta_T_exhaust = Calculate_Delta_T_Exhaust(); +// Propeller calculations. +void FGNewEngine::Do_Prop_Calcs() +{ + float Gear_Ratio = 1.0; + float blade_length; // meters + float forward_velocity; // m/s + float prop_power_consumed_SI; // Watts + float prop_power_consumed_HP; // HP + double J; // advance ratio - dimensionless + double Cp_20; // coefficient of power for 20 degree blade angle + double Cp_25; // coefficient of power for 25 degree blade angle + double Cp; // Our actual coefficient of power + double neta_prop_20; + double neta_prop_25; + double neta_prop; // prop efficiency - // cout << "T_amb " << T_amb; - // cout << " dT exhaust = " << delta_T_exhaust; - - EGT = T_amb + delta_T_exhaust; - - //The above gives the exhaust temperature immediately prior to leaving the combustion chamber - //Now derate to give a more realistic figure as measured downstream - //For now we will aim for a peak of around 400 degC (750 degF) - - EGT *= 0.444 + ((0.544 - 0.444) * Percentage_Power / 100.0); - - EGT_degF = (EGT * 1.8) - 459.67; - - //cout << " EGT = " << EGT << " degK " << EGT_degF << " degF";// << '\n'; - -//*************************************************************************************** -// Calculate Cylinder Head Temperature - -/* DCL 27/10/00 - -This is a somewhat rough first attempt at modelling cylinder head temperature. The cylinder head -is assumed to be at uniform temperature. Obviously this is incorrect, but it simplifies things a -lot, and we're just looking for the behaviour of CHT to be correct. Energy transfer to the cylinder -head is assumed to be one third of the energy released by combustion at all conditions. This is a -reasonable estimate, although obviously in real life it varies with different conditions and possibly -with CHT itself. I've split energy transfer from the cylinder head into 2 terms - free convection - -ie convection to stationary air, and forced convection, ie convection into flowing air. The basic -free convection equation is: dqdt = -hAdT Since we don't know A and are going to set h quite arbitarily -anyway I've knocked A out and just wrapped it up in h - the only real significance is that the units -of h will be different but that dosn't really matter to us anyway. In addition, we have the problem -that the prop model I'm currently using dosn't model the backwash from the prop which will add to the -velocity of the cooling air when the prop is turning, so I've added an extra term to try and cope -with this. - -In real life, forced convection equations are genarally empirically derived, and are quite complicated -and generally contain such things as the Reynolds and Nusselt numbers to various powers. The best -course of action would probably to find an empirical correlation from the literature for a similar -situation and try and get it to fit well. However, for now I am using my own made up very simple -correlation for the energy transfer from the cylinder head: - -dqdt = -(h1.dT) -(h2.m_dot.dT) -(h3.rpm.dT) - -where dT is the temperature different between the cylinder head and the surrounding air, m_dot is the -mass flow rate of cooling air through an arbitary volume, rpm is the engine speed in rpm (this is the -backwash term), and h1, h2, h3 are co-efficients which we can play with to attempt to get the CHT -behaviour to match real life. - -In order to change the values of CHT that the engine settles down at at various conditions, -have a play with h1, h2 and h3. In order to change the rate of heating/cooling without affecting -equilibrium values alter the cylinder head mass, which is really quite arbitary. Bear in mind that -altering h1, h2 and h3 will also alter the rate of heating or cooling as well as equilibrium values, -but altering the cylinder head mass will only alter the rate. It would I suppose be better to read -the values from file to avoid the necessity for re-compilation every time I change them. - -*/ - //CHT = Calc_CHT( Fuel_Flow, Mixture, IAS); - // cout << "Cylinder Head Temp (F) = " << CHT << endl; - float h1 = -95.0; //co-efficient for free convection - float h2 = -3.95; //co-efficient for forced convection - float h3 = -0.05; //co-efficient for forced convection due to prop backwash - float v_apparent; //air velocity over cylinder head in m/s - float v_dot_cooling_air; - float m_dot_cooling_air; - float temperature_difference; - float arbitary_area = 1.0; - float dqdt_from_combustion; - float dqdt_forced; //Rate of energy transfer to/from cylinder head due to forced convection (Joules) (sign convention: to cylinder head is +ve) - float dqdt_free; //Rate of energy transfer to/from cylinder head due to free convection (Joules) (sign convention: to cylinder head is +ve) - float dqdt_cylinder_head; //Overall energy change in cylinder head - float CpCylinderHead = 800.0; //FIXME - this is a guess - I need to look up the correct value - float MassCylinderHead = 8.0; //Kg - this is a guess - it dosn't have to be absolutely accurate but can be varied to alter the warm-up rate - float HeatCapacityCylinderHead; - float dCHTdt; - - temperature_difference = CHT - T_amb; - - v_apparent = IAS * 0.5144444; //convert from knots to m/s - v_dot_cooling_air = arbitary_area * v_apparent; - m_dot_cooling_air = v_dot_cooling_air * rho_air; - - //Calculate rate of energy transfer to cylinder head from combustion - dqdt_from_combustion = m_dot_fuel * calorific_value_fuel * combustion_efficiency * 0.33; - - //Calculate rate of energy transfer from cylinder head due to cooling NOTE is calculated as rate to but negative - dqdt_forced = (h2 * m_dot_cooling_air * temperature_difference) + (h3 * RPM * temperature_difference); - dqdt_free = h1 * temperature_difference; - - //Calculate net rate of energy transfer to or from cylinder head - dqdt_cylinder_head = dqdt_from_combustion + dqdt_forced + dqdt_free; - - HeatCapacityCylinderHead = CpCylinderHead * MassCylinderHead; - - dCHTdt = dqdt_cylinder_head / HeatCapacityCylinderHead; - - CHT += (dCHTdt * time_step); - - CHT_degF = (CHT * 1.8) - 459.67; - - //cout << " CHT = " << CHT_degF << " degF\n"; - - - // End calculate Cylinder Head Temperature - - -//*************************************************************************************** -// Oil pressure calculation - - // Calculate Oil Pressure - Oil_Pressure = Calc_Oil_Press( Oil_Temp, RPM ); - // cout << "Oil Pressure (PSI) = " << Oil_Pressure << endl; - -//************************************************************************************** -// Now do the Propeller Calculations - - Gear_Ratio = 1.0; - FGProp1_RPS = RPM * Gear_Ratio / 60.0; // Borrow this variable from Phils model for now !! + FGProp1_RPS = RPM * Gear_Ratio / 60.0; angular_velocity_SI = 2.0 * LS_PI * RPM / 60.0; forward_velocity = IAS * 0.514444444444; // Convert to m/s - //cout << "Gear_Ratio = " << Gear_Ratio << '\n'; - //cout << "IAS = " << IAS << '\n'; - //cout << "forward_velocity = " << forward_velocity << '\n'; - //cout << "FGProp1_RPS = " << FGProp1_RPS << '\n'; - //cout << "prop_diameter = " << prop_diameter << '\n'; if(FGProp1_RPS == 0) J = 0; else @@ -782,26 +705,4 @@ the values from file to avoid the necessity for re-compilation every time I chan else prop_thrust = neta_prop * prop_power_consumed_SI / forward_velocity; //TODO - rename forward_velocity to IAS_SI //cout << "prop_thrust = " << prop_thrust << '\n'; - -//****************************************************************************** -// Now do the engine - prop torque balance to calculate final RPM - - //Calculate new RPM from torque balance and inertia. - Torque_Imbalance = Torque_SI - prop_torque; //This gives a +ve value when the engine torque exeeds the prop torque - // (Engine torque is +ve when it acts in the direction of engine revolution, prop torque is +ve when it opposes the direction of engine revolution) - - angular_acceleration = Torque_Imbalance / (engine_inertia + prop_inertia); - angular_velocity_SI += (angular_acceleration * time_step); - // Don't let the engine go into reverse - if(angular_velocity_SI < 0) - angular_velocity_SI = 0; - RPM = (angular_velocity_SI * 60) / (2.0 * LS_PI); - -// if(RPM < 0) -// RPM = 0; - - //DCL - stall the engine if RPM drops below 500 - this is possible if mixture lever is pulled right out -// if(RPM < 500) -// RPM = 0; - } diff --git a/src/FDM/IO360.hxx b/src/FDM/IO360.hxx index 76bd14857..2d6c18cb1 100644 --- a/src/FDM/IO360.hxx +++ b/src/FDM/IO360.hxx @@ -1,10 +1,9 @@ -// Module: 10520c.c -// Author: Phil Schubert -// Date started: 12/03/99 -// Purpose: Models a Continental IO-520-M Engine -// Called by: FGSimExec -// -// Copyright (C) 1999 Philip L. Schubert (philings@ozemail.com.au) +// IO360.hxx - a piston engine model currently for the IO360 engine fitted to the C172 +// but with the potential to model other naturally aspirated piston engines +// given appropriate config input. +// +// Written by David Luff, started 2000. +// Based on code by Phil Schubert, started 1999. // // This program is free software; you can redistribute it and/or // modify it under the terms of the GNU General Public License as @@ -18,34 +17,8 @@ // // You should have received a copy of the GNU General Public License // along with this program; if not, write to the Free Software -// Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA -// 02111-1307, USA. -// -// Further information about the GNU General Public License can also -// be found on the world wide web at http://www.gnu.org. -// -// FUNCTIONAL DESCRIPTION -// ------------------------------------------------------------------------ -// Models a Continental IO-520-M engine. This engine is used in Cessna -// 210, 310, Beechcraft Bonaza and Baron C55. The equations used below -// were determined by a first and second order curve fits using Excel. -// The data is from the Cessna Aircraft Corporations Engine and Flight -// Computer for C310. Part Number D3500-13 -// -// ARGUMENTS -// ------------------------------------------------------------------------ -// -// -// HISTORY -// ------------------------------------------------------------------------ -// 12/03/99 PLS Created -// 07/03/99 PLS Added Calculation of Density, and Prop_Torque -// 07/03/99 PLS Restructered Variables to allow easier implementation -// of Classes -// 15/03/99 PLS Added Oil Pressure, Oil Temperature and CH Temp -// ------------------------------------------------------------------------ -// INCLUDES -// ------------------------------------------------------------------------ +// Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. + #ifndef _IO360_HXX_ #define _IO360_HXX_ @@ -65,109 +38,107 @@ class FGNewEngine { private: - float CONVERT_HP_TO_WATTS; + // These constants should probably be moved eventually float CONVERT_CUBIC_INCHES_TO_METERS_CUBED; + float CONVERT_HP_TO_WATTS; - // Control and environment inputs + // Properties of working fluids + float Cp_air; // J/KgK + float Cp_fuel; // J/KgK + float calorific_value_fuel; // W/Kg + float rho_fuel; // kg/m^3 + float rho_air; // kg/m^3 + + // environment inputs + float p_amb_sea_level; // Sea level ambient pressure in Pascals + float p_amb; // Ambient pressure at current altitude in Pascals + float T_amb; // ditto deg Kelvin + + // Control inputs + float Throttle_Lever_Pos; // 0 = Closed, 100 = Fully Open + float Propeller_Lever_Pos; // 0 = Full Course 100 = Full Fine + float Mixture_Lever_Pos; // 0 = Idle Cut Off 100 = Full Rich + + //misc float IAS; - // 0 = Closed, 100 = Fully Open - float Throttle_Lever_Pos; - // 0 = Full Course 100 = Full Fine - float Propeller_Lever_Pos; - // 0 = Idle Cut Off 100 = Full Rich - float Mixture_Lever_Pos; + double time_step; - // Engine Specific Variables used by this program that have limits. - // Will be set in a parameter file to be read in to create - // and instance for each engine. - float Max_Manifold_Pressure; //will be lower than ambient pressure for a non turbo/super charged engine due to losses through the throttle. This is the sea level full throttle value. - float Min_Manifold_Pressure; //Closed throttle valueat idle - governed by the idle bypass valve - float Max_RPM; - float Min_RPM; - float Max_Fuel_Flow; - float Mag_Derate_Percent; - float MaxHP; - float Gear_Ratio; + // Engine Specific Variables that should be read in from a config file + float MaxHP; // Horsepower + float displacement; // Cubic inches + float displacement_SI; //m^3 (derived from above rather than read in) + float engine_inertia; //kg.m^2 + float prop_inertia; //kg.m^2 + float Max_Fuel_Flow; // Units??? Do we need this variable any more?? - // Initialise Engine Variables used by this instance + // Engine specific variables that maybe should be read in from config but are pretty generic and won't vary much for a naturally aspirated piston engine. + float Max_Manifold_Pressure; // inches Hg - typical manifold pressure at full throttle and typical max rpm + float Min_Manifold_Pressure; // inches Hg - typical manifold pressure at idle (closed throttle) + float Max_RPM; // rpm - this is really a bit of a hack and could be make redundant if the prop model works properly and takes tips at speed of sound into account. + float Min_RPM; // rpm - possibly redundant ??? + float Mag_Derate_Percent; // Percentage reduction in power when mags are switched from 'both' to either 'L' or 'R' + float Gear_Ratio; // Gearing between engine and propellor + float n_R; // Number of cycles per power stroke - 2 for a 4 stroke engine. + + // Engine Variables not read in from config file + float RPM; // rpm float Percentage_Power; // Power output as percentage of maximum power output - float Manifold_Pressure; // Inches - float RPM; - float Fuel_Flow_gals_hr; // gals/hour - float Torque; - float CHT; // Cylinder head temperature deg K + float Manifold_Pressure; // Inches Hg + float Fuel_Flow_gals_hr; // USgals/hour + float Torque_lbft; // lb-ft + float Torque_SI; // Nm + float CHT_degK; // Cylinder head temperature deg K float CHT_degF; // Ditto in deg Fahrenheit - float EGT; // Exhaust gas temperature deg K - float EGT_degF; // Exhaust gas temperature deg Fahrenheit float Mixture; float Oil_Pressure; // PSI float Oil_Temp; // Deg C + float current_oil_temp; // deg K + /**** one of these is superfluous !!!!***/ float HP; // Current power output in HP float Power_SI; // Current power output in Watts - float Torque_SI; // Torque in Nm - float Torque_FMEP; // The component of Engine torque due to FMEP (Nm) float RPS; - float Torque_Imbalance; - bool running; //flag to indicate the engine is running self sustaining - bool cranking; //flag to indicate the engine is being cranked - bool spark; //flag to indicate a spark is available - bool fuel; //flag to indicate fuel is available - - //DCL + float angular_velocity_SI; // rad/s + float Torque_FMEP; // The component of Engine torque due to FMEP (Nm) + float Torque_Imbalance; // difference between engine and prop torque + float EGT; // Exhaust gas temperature deg K + float EGT_degF; // Exhaust gas temperature deg Fahrenheit float volumetric_efficiency; float combustion_efficiency; - float equivalence_ratio; - float v_dot_air; - float m_dot_air; - float m_dot_fuel; - float swept_volume; - float True_Manifold_Pressure; //in Hg - float rho_air_manifold; - float R_air; - float p_amb_sea_level; // Pascals - float p_amb; // Pascals - float T_amb; // deg Kelvin - float calorific_value_fuel; - float rho_air; - float rho_fuel; // kg/m^3 - float thi_sea_level; - float delta_T_exhaust; - float displacement; // Engine displacement in cubic inches - to be read in from config file for each engine - float displacement_SI; // ditto in meters cubed - float Cp_air; // J/KgK - float Cp_fuel; // J/KgK - float heat_capacity_exhaust; - float enthalpy_exhaust; - float Percentage_of_best_power_mixture_power; + float equivalence_ratio; // ratio of stoichiometric AFR over actual AFR + float thi_sea_level; // the equivalence ratio we would be running at assuming sea level air denisity in the manifold + float v_dot_air; // volume flow rate of air into engine - m^3/s + float m_dot_air; // mass flow rate of air into engine - kg/s + float m_dot_fuel; // mass flow rate of fuel into engine - kg/s + float swept_volume; // total engine swept volume - m^3 + /********* swept volume or the geometry used to calculate it should be in the config read section surely ??? ******/ + float True_Manifold_Pressure; //in Hg - actual manifold gauge pressure + float rho_air_manifold; // denisty of air in the manifold - kg/m^3 + float R_air; // Gas constant of air (287) UNITS??? + float delta_T_exhaust; // Temperature change of exhaust this time step - degK + float heat_capacity_exhaust; // exhaust gas specific heat capacity - J/kgK + float enthalpy_exhaust; // Enthalpy at current exhaust gas conditions - UNITS??? + float Percentage_of_best_power_mixture_power; // Current power as a percentage of what power we would have at the same conditions but at best power mixture float abstract_mixture; //temporary hack - float engine_inertia; //kg.m^2 - float prop_inertia; //kg.m^2 float angular_acceleration; //rad/s^2 - float n_R; //Number of cycles per power stroke float FMEP; //Friction Mean Effective Pressure (Pa) - double time_step; + + // Various bits of housekeeping describing the engines state. + bool running; // flag to indicate the engine is running self sustaining + bool cranking; // flag to indicate the engine is being cranked + int crank_counter; // Number of iterations that the engine has been cranked non-stop + bool spark; // flag to indicate a spark is available + bool fuel; // flag to indicate fuel is available // Propellor Variables - float FGProp1_Thrust; - float FGProp1_RPS; - float FGProp1_Blade_Angle; + float FGProp1_RPS; // rps float prop_torque; // Nm float prop_thrust; // Newtons - float blade_length; // meters - float forward_velocity; // m/s - float angular_velocity_SI; // rad/s - float prop_power_consumed_SI; // Watts - float prop_power_consumed_HP; // HP - double prop_diameter; // meters - double J; // advance ratio - dimensionless - double Cp_20; // coefficient of power for 20 degree blade angle - double Cp_25; // coefficient of power for 25 degree blade angle - double Cp; // Our actual coefficient of power - double blade_angle; // degrees - double neta_prop_20; - double neta_prop_25; - double neta_prop; // prop efficiency + double prop_diameter; // meters + double blade_angle; // degrees + +// MEMBER FUNCTIONS + // Calculate Engine RPM based on Propellor Lever Position float Calc_Engine_RPM(float Position); @@ -181,18 +152,32 @@ private: // Calculate percentage of best power mixture power based on equivalence ratio float Power_Mixture_Correlation(float thi_actual); - // Calculate exhaust gas temperature rise + // Calculate exhaust gas temperature change float Calculate_Delta_T_Exhaust(void); + // Calculate cylinder head temperature + float FGNewEngine::Calc_CHT(float CHT); + + void FGNewEngine::Calc_EGT(void); + + // Calculate fuel flow in gals/hr + void FGNewEngine::Calc_Fuel_Flow_Gals_Hr(void); + + // Calculate current percentage power + void FGNewEngine::Calc_Percentage_Power(bool mag_left, bool mag_right); + // Calculate Oil Temperature - float Calc_Oil_Temp (float Fuel_Flow, float Mixture, float IAS); + float Calc_Oil_Temp (float oil_temp); // Calculate Oil Pressure float Calc_Oil_Press (float Oil_Temp, float Engine_RPM); + // Propeller calculations. + void FGNewEngine::Do_Prop_Calcs(void); + public: - ofstream outfile; +// ofstream outfile; //constructor FGNewEngine() { @@ -234,9 +219,6 @@ public: // accessors inline float get_RPM() const { return RPM; } inline float get_Manifold_Pressure() const { return True_Manifold_Pressure; } - inline float get_FGProp1_Thrust() const { return FGProp1_Thrust; } - inline float get_FGProp1_Blade_Angle() const { return FGProp1_Blade_Angle; } - // inline float get_Rho() const { return Rho; } inline float get_MaxHP() const { return MaxHP; } inline float get_Percentage_Power() const { return Percentage_Power; } @@ -245,6 +227,7 @@ public: inline float get_prop_thrust_SI() const { return prop_thrust; } inline float get_prop_thrust_lbs() const { return (prop_thrust * 0.2248); } inline float get_fuel_flow_gals_hr() const { return (Fuel_Flow_gals_hr); } + inline float get_oil_temp() const { return ((current_oil_temp * 1.8) - 459.67); } };