From b7ebc7d78dbf320d432bad3ed0e384b2e1308eed Mon Sep 17 00:00:00 2001 From: ehofman Date: Sun, 30 Nov 2008 10:44:29 +0000 Subject: [PATCH] Sync. with JSBSim CVS --- src/FDM/JSBSim/FGState.cpp | 2 +- src/FDM/JSBSim/JSBSim.cxx | 4 + src/FDM/JSBSim/input_output/FGXMLElement.cpp | 12 +- src/FDM/JSBSim/input_output/FGXMLElement.h | 6 + src/FDM/JSBSim/models/FGAtmosphere.cpp | 101 ++++++--- src/FDM/JSBSim/models/FGAtmosphere.h | 127 +++++++---- src/FDM/JSBSim/models/FGAuxiliary.cpp | 8 +- src/FDM/JSBSim/models/FGFCS.cpp | 8 +- src/FDM/JSBSim/models/FGMassBalance.cpp | 16 +- src/FDM/JSBSim/models/FGMassBalance.h | 39 ++-- src/FDM/JSBSim/models/FGOutput.cpp | 4 +- src/FDM/JSBSim/models/FGPropulsion.cpp | 6 +- src/FDM/JSBSim/models/atmosphere/FGMSIS.cpp | 6 +- src/FDM/JSBSim/models/atmosphere/FGMars.cpp | 193 +---------------- src/FDM/JSBSim/models/atmosphere/FGMars.h | 103 +-------- src/FDM/JSBSim/models/propulsion/FGElectric.h | 3 +- src/FDM/JSBSim/models/propulsion/FGEngine.cpp | 44 ++-- src/FDM/JSBSim/models/propulsion/FGEngine.h | 40 ++-- src/FDM/JSBSim/models/propulsion/FGNozzle.cpp | 47 +--- src/FDM/JSBSim/models/propulsion/FGNozzle.h | 19 +- src/FDM/JSBSim/models/propulsion/FGPiston.cpp | 69 +++--- src/FDM/JSBSim/models/propulsion/FGPiston.h | 13 +- src/FDM/JSBSim/models/propulsion/FGRocket.cpp | 202 ++++++++++++++---- src/FDM/JSBSim/models/propulsion/FGRocket.h | 69 ++++-- src/FDM/JSBSim/models/propulsion/FGTank.cpp | 73 ++++++- src/FDM/JSBSim/models/propulsion/FGTank.h | 22 +- .../JSBSim/models/propulsion/FGTurbine.cpp | 7 +- src/FDM/JSBSim/models/propulsion/FGTurbine.h | 2 +- .../JSBSim/models/propulsion/FGTurboProp.cpp | 5 +- 29 files changed, 640 insertions(+), 610 deletions(-) diff --git a/src/FDM/JSBSim/FGState.cpp b/src/FDM/JSBSim/FGState.cpp index f5fc784f2..cc6bbfa57 100644 --- a/src/FDM/JSBSim/FGState.cpp +++ b/src/FDM/JSBSim/FGState.cpp @@ -96,7 +96,7 @@ void FGState::Initialize(FGInitialCondition *FGIC) FGIC->GetWindDFpsIC() ); FGColumnVector3 vAeroUVW; - vAeroUVW = Propagate->GetUVW() + Propagate->GetTl2b()*Atmosphere->GetWindNED(); + vAeroUVW = Propagate->GetUVW() + Propagate->GetTl2b()*Atmosphere->GetTotalWindNED(); double alpha, beta; if (vAeroUVW(eW) != 0.0) diff --git a/src/FDM/JSBSim/JSBSim.cxx b/src/FDM/JSBSim/JSBSim.cxx index 9ae847043..f63e106ad 100644 --- a/src/FDM/JSBSim/JSBSim.cxx +++ b/src/FDM/JSBSim/JSBSim.cxx @@ -750,6 +750,10 @@ bool FGJSBsim::copy_from_JSBsim() node->setDoubleValue("oil-temperature-degf", eng->getOilTemp_degF()); node->setDoubleValue("oil-pressure-psi", eng->getOilPressure_psi()); node->setDoubleValue("mp-osi", eng->getManifoldPressure_inHg()); + // NOTE: mp-osi is not in ounces per square inch. + // This error is left for reasons of backwards compatibility with + // existing FlightGear sound and instrument configurations. + node->setDoubleValue("mp-inhg", eng->getManifoldPressure_inHg()); node->setDoubleValue("cht-degf", eng->getCylinderHeadTemp_degF()); node->setDoubleValue("rpm", eng->getRPM()); } // end FGPiston code block diff --git a/src/FDM/JSBSim/input_output/FGXMLElement.cpp b/src/FDM/JSBSim/input_output/FGXMLElement.cpp index 697577f5b..74983440e 100755 --- a/src/FDM/JSBSim/input_output/FGXMLElement.cpp +++ b/src/FDM/JSBSim/input_output/FGXMLElement.cpp @@ -65,6 +65,8 @@ Element::Element(string nm) convert["FT2"]["M2"] = 1.0/convert["M2"]["FT2"]; convert["M2"]["IN2"] = convert["M"]["IN"]*convert["M"]["IN"]; convert["IN2"]["M2"] = 1.0/convert["M2"]["IN2"]; + convert["FT2"]["IN2"] = 144.0; + convert["IN2"]["FT2"] = 1.0/convert["FT2"]["IN2"]; // Volume convert["IN3"]["CC"] = 16.387064; convert["CC"]["IN3"] = 1.0/convert["IN3"]["CC"]; @@ -128,6 +130,9 @@ Element::Element(string nm) convert["PA"]["LBS/FT2"] = 1.0/convert["LBS/FT2"]["PA"]; // Mass flow convert["KG/MIN"]["LBS/MIN"] = convert["KG"]["LBS"]; + // Fuel Consumption + convert["LBS/HP*HR"]["KG/KW*HR"] = 0.6083; + convert["KG/KW*HR"]["LBS/HP*HR"] = 1.0/convert["LBS/HP*HR"]["KG/KW*HR"]; // Length convert["M"]["M"] = 1.00; @@ -187,6 +192,9 @@ Element::Element(string nm) convert["LBS/SEC"]["LBS/SEC"] = 1.00; convert["KG/MIN"]["KG/MIN"] = 1.0; convert["LBS/MIN"]["LBS/MIN"] = 1.0; + // Fuel Consumption + convert["LBS/HP*HR"]["LBS/HP*HR"] = 1.0; + convert["KG/KW*HR"]["KG/KW*HR"] = 1.0; } //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% @@ -389,12 +397,12 @@ double Element::FindElementValueAsNumberConvertFromTo( string el, string target_units) { Element* element = FindElement(el); - + if (!element) { cerr << "Attempting to get non-existent element " << el << endl; exit(0); } - + if (!supplied_units.empty()) { if (convert.find(supplied_units) == convert.end()) { cerr << endl << "Supplied unit: \"" << supplied_units << "\" does not exist (typo?). Add new unit" diff --git a/src/FDM/JSBSim/input_output/FGXMLElement.h b/src/FDM/JSBSim/input_output/FGXMLElement.h index acfb07436..90bfa58f7 100755 --- a/src/FDM/JSBSim/input_output/FGXMLElement.h +++ b/src/FDM/JSBSim/input_output/FGXMLElement.h @@ -90,6 +90,8 @@ CLASS DOCUMENTATION - convert["LBS"]["N"] = 1.0/convert["N"]["LBS"]; - convert["KTS"]["FT/SEC"] = ktstofps; - convert["KG/MIN"]["LBS/MIN"] = convert["KG"]["LBS"]; + - convert["LBS/HP*HR"]["KG/KW*HR"] = 0.6083; + - convert["KG/KW*HR"]["LBS/HP*HR"] = 1/convert["LBS/HP*HR"]["KG/KW*HR"]; - convert["M"]["M"] = 1.00; - convert["FT"]["FT"] = 1.00; @@ -115,6 +117,8 @@ CLASS DOCUMENTATION - convert["FT/SEC"]["FT/SEC"] = 1.0; - convert["KG/MIN"]["KG/MIN"] = 1.0; - convert["LBS/MIN"]["LBS/MIN"] = 1.0; + - convert["LBS/HP*HR"]["LBS/HP*HR"] = 1.0; + - convert["KG/KW*HR"]["KG/KW*HR"] = 1.0; Where: - N = newtons @@ -130,6 +134,8 @@ CLASS DOCUMENTATION - DEG = degrees - RAD = radians - WATTS = watts + - HP = horsepower + - HR = hour @author Jon S. Berndt @version $Id$ diff --git a/src/FDM/JSBSim/models/FGAtmosphere.cpp b/src/FDM/JSBSim/models/FGAtmosphere.cpp index dda86fc13..88ffbd5e0 100644 --- a/src/FDM/JSBSim/models/FGAtmosphere.cpp +++ b/src/FDM/JSBSim/models/FGAtmosphere.cpp @@ -90,8 +90,8 @@ FGAtmosphere::FGAtmosphere(FGFDMExec* fdmex) : FGModel(fdmex) T_dev_sl = T_dev = delta_T = 0.0; StandardTempOnly = false; first_pass = true; - vGustNED(1) = vGustNED(2) = vGustNED(3) = 0.0; bgustSet = false; - vTurbulence(1) = vTurbulence(2) = vTurbulence(3) = 0.0; + vGustNED.InitMatrix(); + vTurbulenceNED.InitMatrix(); bind(); Debug(0); @@ -256,17 +256,18 @@ void FGAtmosphere::Calculate(double altitude) //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% // Calculate parameters derived from T, P and rho +// Sum gust and turbulence values in NED frame into the wind vector. void FGAtmosphere::CalculateDerived(void) { T_dev = (*temperature) - GetTemperature(h); density_altitude = h + T_dev * 66.7; - if (turbType == ttStandard || ttCulp) { - Turbulence(); - vWindNED += vGustNED + vTurbulence; - } - if (vWindNED(1) != 0.0) psiw = atan2( vWindNED(2), vWindNED(1) ); + if (turbType == ttStandard || ttCulp) Turbulence(); + + vTotalWindNED = vWindNED + vGustNED + vTurbulenceNED; + + if (vWindNED(eX) != 0.0) psiw = atan2( vWindNED(eY), vWindNED(eX) ); if (psiw < 0) psiw += 2*M_PI; soundspeed = sqrt(SHRatio*Reng*(*temperature)); @@ -326,6 +327,36 @@ static inline double square_signed (double value) //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +void FGAtmosphere::SetWindspeed(double speed) +{ + if (vWindNED.Magnitude() == 0.0) { + psiw = 0.0; + vWindNED(eNorth) = speed; + } else { + vWindNED(eNorth) = speed * cos(psiw); + vWindNED(eEast) = speed * sin(psiw); + vWindNED(eDown) = 0.0; + } +} + +//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +double FGAtmosphere::GetWindspeed(void) const +{ + return vWindNED.Magnitude(); +} + +//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +void FGAtmosphere::SetWindPsi(double dir) +{ + double mag = GetWindspeed(); + psiw = dir; + SetWindspeed(mag); +} + +//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + void FGAtmosphere::Turbulence(void) { switch (turbType) { @@ -359,10 +390,10 @@ void FGAtmosphere::Turbulence(void) // Diminish turbulence within three wingspans // of the ground - vTurbulence = TurbGain * Magnitude * vDirection; + vTurbulenceNED = TurbGain * Magnitude * vDirection; double HOverBMAC = Auxiliary->GetHOverBMAC(); if (HOverBMAC < 3.0) - vTurbulence *= (HOverBMAC / 3.0) * (HOverBMAC / 3.0); + vTurbulenceNED *= (HOverBMAC / 3.0) * (HOverBMAC / 3.0); // I don't believe these next two statements calculate the proper gradient over // the aircraft body. One reason is because this has no relationship with the @@ -395,8 +426,8 @@ void FGAtmosphere::Turbulence(void) // that we've used them to calculate // moments. // Why? (JSB) -// vTurbulence(eX) = 0.0; -// vTurbulence(eY) = 0.0; +// vTurbulenceNED(eX) = 0.0; +// vTurbulenceNED(eY) = 0.0; break; } @@ -426,7 +457,7 @@ void FGAtmosphere::Turbulence(void) vDirection.Normalize(); - vTurbulence = TurbGain*Magnitude * vDirection; + vTurbulenceNED = TurbGain*Magnitude * vDirection; vTurbulenceGrad = TurbGain*MagnitudeAccel * vDirection; vBodyTurbGrad = Propagate->GetTl2b()*vTurbulenceGrad; @@ -473,19 +504,19 @@ void FGAtmosphere::Turbulence(void) // max vertical wind speed in fps, corresponds to TurbGain = 1.0 double max_vs = 40; - vTurbulence(1) = vTurbulence(2) = vTurbulence(3) = 0.0; + vTurbulenceNED(1) = vTurbulenceNED(2) = vTurbulenceNED(3) = 0.0; double delta = strength * max_vs * TurbGain * (1-Rhythmicity) * spike; // Vertical component of turbulence. - vTurbulence(3) = sinewave * max_vs * TurbGain * Rhythmicity; - vTurbulence(3)+= delta; + vTurbulenceNED(3) = sinewave * max_vs * TurbGain * Rhythmicity; + vTurbulenceNED(3)+= delta; double HOverBMAC = Auxiliary->GetHOverBMAC(); if (HOverBMAC < 3.0) - vTurbulence(3) *= HOverBMAC * 0.3333; + vTurbulenceNED(3) *= HOverBMAC * 0.3333; // Yaw component of turbulence. - vTurbulence(1) = sin( delta * 3.0 ); - vTurbulence(2) = cos( delta * 3.0 ); + vTurbulenceNED(1) = sin( delta * 3.0 ); + vTurbulenceNED(2) = cos( delta * 3.0 ); // Roll component of turbulence. Clockwise vortex causes left roll. vTurbPQR(eP) += delta * 0.04; @@ -541,6 +572,28 @@ void FGAtmosphere::bind(void) PropertyManager->Tie("atmosphere/delta-T", this, &FGAtmosphere::GetDeltaT, &FGAtmosphere::SetDeltaT); PropertyManager->Tie("atmosphere/T-sl-dev-F", this, &FGAtmosphere::GetSLTempDev, &FGAtmosphere::SetSLTempDev); PropertyManager->Tie("atmosphere/density-altitude", this, &FGAtmosphere::GetDensityAltitude); + + PropertyManager->Tie("atmosphere/wind-north-fps", this, eNorth, (PMF)&FGAtmosphere::GetWindNED, + (PMFd)&FGAtmosphere::SetWindNED); + PropertyManager->Tie("atmosphere/wind-east-fps", this, eEast, (PMF)&FGAtmosphere::GetWindNED, + (PMFd)&FGAtmosphere::SetWindNED); + PropertyManager->Tie("atmosphere/wind-down-fps", this, eDown, (PMF)&FGAtmosphere::GetWindNED, + (PMFd)&FGAtmosphere::SetWindNED); + PropertyManager->Tie("atmosphere/wind-from-cw", this, &FGAtmosphere::GetWindFromClockwise, + &FGAtmosphere::SetWindFromClockwise); + PropertyManager->Tie("atmosphere/wind-mag-fps", this, &FGAtmosphere::GetWindspeed, + &FGAtmosphere::SetWindspeed); + PropertyManager->Tie("atmosphere/total-wind-north-fps", this, eNorth, (PMF)&FGAtmosphere::GetTotalWindNED); + PropertyManager->Tie("atmosphere/total-wind-east-fps", this, eEast, (PMF)&FGAtmosphere::GetTotalWindNED); + PropertyManager->Tie("atmosphere/total-wind-down-fps", this, eDown, (PMF)&FGAtmosphere::GetTotalWindNED); + + PropertyManager->Tie("atmosphere/gust-north-fps", this, eNorth, (PMF)&FGAtmosphere::GetGustNED, + (PMFd)&FGAtmosphere::SetGustNED); + PropertyManager->Tie("atmosphere/gust-east-fps", this, eEast, (PMF)&FGAtmosphere::GetGustNED, + (PMFd)&FGAtmosphere::SetGustNED); + PropertyManager->Tie("atmosphere/gust-down-fps", this, eDown, (PMF)&FGAtmosphere::GetGustNED, + (PMFd)&FGAtmosphere::SetGustNED); + PropertyManager->Tie("atmosphere/p-turb-rad_sec", this,1, (PMF)&FGAtmosphere::GetTurbPQR); PropertyManager->Tie("atmosphere/q-turb-rad_sec", this,2, (PMF)&FGAtmosphere::GetTurbPQR); PropertyManager->Tie("atmosphere/r-turb-rad_sec", this,3, (PMF)&FGAtmosphere::GetTurbPQR); @@ -548,14 +601,6 @@ void FGAtmosphere::bind(void) PropertyManager->Tie("atmosphere/turb-gain", this, &FGAtmosphere::GetTurbGain, &FGAtmosphere::SetTurbGain); PropertyManager->Tie("atmosphere/turb-rhythmicity", this, &FGAtmosphere::GetRhythmicity, &FGAtmosphere::SetRhythmicity); - PropertyManager->Tie("atmosphere/gust-north-fps", this,1, (PMF)&FGAtmosphere::GetGustNED, - (PMFd)&FGAtmosphere::SetGustNED); - PropertyManager->Tie("atmosphere/gust-east-fps", this,2, (PMF)&FGAtmosphere::GetGustNED, - (PMFd)&FGAtmosphere::SetGustNED); - PropertyManager->Tie("atmosphere/gust-down-fps", this,3, (PMF)&FGAtmosphere::GetGustNED, - (PMFd)&FGAtmosphere::SetGustNED); - PropertyManager->Tie("atmosphere/wind-from-cw", this, &FGAtmosphere::GetWindFromClockwise, - &FGAtmosphere::SetWindFromClockwise); } //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% @@ -598,14 +643,14 @@ void FGAtmosphere::Debug(int from) if (debug_lvl & 128) { // Turbulence if (first_pass && from == 2) { first_pass = false; - cout << "vTurbulence(X), vTurbulence(Y), vTurbulence(Z), " + cout << "vTurbulenceNED(X), vTurbulenceNED(Y), vTurbulenceNED(Z), " << "vTurbulenceGrad(X), vTurbulenceGrad(Y), vTurbulenceGrad(Z), " << "vDirection(X), vDirection(Y), vDirection(Z), " << "Magnitude, " << "vTurbPQR(P), vTurbPQR(Q), vTurbPQR(R), " << endl; } if (from == 2) { - cout << vTurbulence << ", " << vTurbulenceGrad << ", " << vDirection << ", " << Magnitude << ", " << vTurbPQR << endl; + cout << vTurbulenceNED << ", " << vTurbulenceGrad << ", " << vDirection << ", " << Magnitude << ", " << vTurbPQR << endl; } } if (debug_lvl & 64) { diff --git a/src/FDM/JSBSim/models/FGAtmosphere.h b/src/FDM/JSBSim/models/FGAtmosphere.h index 4434c3179..ed7bfcf93 100644 --- a/src/FDM/JSBSim/models/FGAtmosphere.h +++ b/src/FDM/JSBSim/models/FGAtmosphere.h @@ -86,10 +86,10 @@ public: enum tType {ttStandard, ttBerndt, ttCulp, ttNone} turbType; /// Returns the temperature in degrees Rankine. - inline double GetTemperature(void) const {return *temperature;} + double GetTemperature(void) const {return *temperature;} /** Returns the density in slugs/ft^3. This function may only be used if Run() is called first. */ - inline double GetDensity(void) const {return *density;} + double GetDensity(void) const {return *density;} /// Returns the pressure in psf. double GetPressure(void) const {return *pressure;} /// Returns the standard pressure at a specified altitude @@ -99,25 +99,25 @@ public: /// Returns the standard density at a specified altitude double GetDensity(double altitude); /// Returns the speed of sound in ft/sec. - inline double GetSoundSpeed(void) const {return soundspeed;} + double GetSoundSpeed(void) const {return soundspeed;} /// Returns the sea level temperature in degrees Rankine. - inline double GetTemperatureSL(void) const { return SLtemperature; } + double GetTemperatureSL(void) const { return SLtemperature; } /// Returns the sea level density in slugs/ft^3 - inline double GetDensitySL(void) const { return SLdensity; } + double GetDensitySL(void) const { return SLdensity; } /// Returns the sea level pressure in psf. - inline double GetPressureSL(void) const { return SLpressure; } + double GetPressureSL(void) const { return SLpressure; } /// Returns the sea level speed of sound in ft/sec. - inline double GetSoundSpeedSL(void) const { return SLsoundspeed; } + double GetSoundSpeedSL(void) const { return SLsoundspeed; } /// Returns the ratio of at-altitude temperature over the sea level value. - inline double GetTemperatureRatio(void) const { return (*temperature)*rSLtemperature; } + double GetTemperatureRatio(void) const { return (*temperature)*rSLtemperature; } /// Returns the ratio of at-altitude density over the sea level value. - inline double GetDensityRatio(void) const { return (*density)*rSLdensity; } + double GetDensityRatio(void) const { return (*density)*rSLdensity; } /// Returns the ratio of at-altitude pressure over the sea level value. - inline double GetPressureRatio(void) const { return (*pressure)*rSLpressure; } + double GetPressureRatio(void) const { return (*pressure)*rSLpressure; } /// Returns the ratio of at-altitude sound speed over the sea level value. - inline double GetSoundSpeedRatio(void) const { return soundspeed*rSLsoundspeed; } + double GetSoundSpeedRatio(void) const { return soundspeed*rSLsoundspeed; } /// Tells the simulator to use an externally calculated atmosphere model. void UseExternal(void); @@ -127,66 +127,100 @@ public: bool External(void) { return useExternal; } /// Provides the external atmosphere model with an interface to set the temperature. - inline void SetExTemperature(double t) { exTemperature=t; } + void SetExTemperature(double t) { exTemperature=t; } /// Provides the external atmosphere model with an interface to set the density. - inline void SetExDensity(double d) { exDensity=d; } + void SetExDensity(double d) { exDensity=d; } /// Provides the external atmosphere model with an interface to set the pressure. - inline void SetExPressure(double p) { exPressure=p; } + void SetExPressure(double p) { exPressure=p; } /// Sets the temperature deviation at sea-level in degrees Fahrenheit - inline void SetSLTempDev(double d) { T_dev_sl = d; } + void SetSLTempDev(double d) { T_dev_sl = d; } /// Gets the temperature deviation at sea-level in degrees Fahrenheit - inline double GetSLTempDev(void) const { return T_dev_sl; } + double GetSLTempDev(void) const { return T_dev_sl; } /// Sets the current delta-T in degrees Fahrenheit - inline void SetDeltaT(double d) { delta_T = d; } + void SetDeltaT(double d) { delta_T = d; } /// Gets the current delta-T in degrees Fahrenheit - inline double GetDeltaT(void) const { return delta_T; } + double GetDeltaT(void) const { return delta_T; } /// Gets the at-altitude temperature deviation in degrees Fahrenheit - inline double GetTempDev(void) const { return T_dev; } + double GetTempDev(void) const { return T_dev; } /// Gets the density altitude in feet - inline double GetDensityAltitude(void) const { return density_altitude; } + double GetDensityAltitude(void) const { return density_altitude; } + + // TOTAL WIND access functions (wind + gust + turbulence) + + /// Retrieves the total wind components in NED frame. + FGColumnVector3& GetTotalWindNED(void) { return vTotalWindNED; } + + /// Retrieves a total wind component in NED frame. + double GetTotalWindNED(int idx) const {return vTotalWindNED(idx);} + + // WIND access functions /// Sets the wind components in NED frame. - inline void SetWindNED(double wN, double wE, double wD) { vWindNED(1)=wN; vWindNED(2)=wE; vWindNED(3)=wD;} + void SetWindNED(double wN, double wE, double wD) { vWindNED(1)=wN; vWindNED(2)=wE; vWindNED(3)=wD;} + + /// Sets a wind component in NED frame. + void SetWindNED(int idx, double wind) { vWindNED(idx)=wind;} /// Retrieves the wind components in NED frame. - inline FGColumnVector3& GetWindNED(void) { return vWindNED; } + FGColumnVector3& GetWindNED(void) { return vWindNED; } - /// Sets gust components in NED frame. - inline void SetGustNED(int idx, double gust) { vGustNED(idx)=gust;} + /// Retrieves a wind component in NED frame. + double GetWindNED(int idx) const {return vWindNED(idx);} - /// Retrieves the gust components in NED frame. - inline double GetGustNED(int idx) const {return vGustNED(idx);} + /** Retrieves the direction that the wind is coming from. + The direction is defined as north=0 and increases counterclockwise. + The wind heading is returned in radians.*/ + double GetWindPsi(void) const { return psiw; } - /// Retrieves the gust components in NED frame. - inline FGColumnVector3& GetGustNED(void) {return vGustNED;} + /** Sets the direction that the wind is coming from. + The direction is defined as north=0 and increases counterclockwise to 2*pi (radians). The + vertical component of wind is assumed to be zero - and is forcibly set to zero. This function + sets the vWindNED vector components based on the supplied direction. The magnitude of + the wind set in the vector is preserved (assuming the vertical component is non-zero). + @param dir wind direction in the horizontal plane, in radians.*/ + void SetWindPsi(double dir); + + void SetWindspeed(double speed); + + double GetWindspeed(void) const; - /** Retrieves the wind direction. The direction is defined as north=0 and - increases counterclockwise. The wind heading is returned in radians.*/ - inline double GetWindPsi(void) const { return psiw; } + // GUST access functions + + /// Sets a gust component in NED frame. + void SetGustNED(int idx, double gust) { vGustNED(idx)=gust;} + + /// Sets the gust components in NED frame. + void SetGustNED(double gN, double gE, double gD) { vGustNED(eNorth)=gN; vGustNED(eEast)=gE; vGustNED(eDown)=gD;} + + /// Retrieves a gust component in NED frame. + double GetGustNED(int idx) const {return vGustNED(idx);} + + /// Retrieves the gust components in NED frame. + FGColumnVector3& GetGustNED(void) {return vGustNED;} /** Turbulence models available: ttStandard, ttBerndt, ttCulp, ttNone */ - inline void SetTurbType(tType tt) {turbType = tt;} - inline tType GetTurbType() const {return turbType;} + void SetTurbType(tType tt) {turbType = tt;} + tType GetTurbType() const {return turbType;} - inline void SetTurbGain(double tg) {TurbGain = tg;} - inline double GetTurbGain() const {return TurbGain;} + void SetTurbGain(double tg) {TurbGain = tg;} + double GetTurbGain() const {return TurbGain;} - inline void SetTurbRate(double tr) {TurbRate = tr;} - inline double GetTurbRate() const {return TurbRate;} + void SetTurbRate(double tr) {TurbRate = tr;} + double GetTurbRate() const {return TurbRate;} - inline void SetRhythmicity(double r) {Rhythmicity=r;} - inline double GetRhythmicity() const {return Rhythmicity;} + void SetRhythmicity(double r) {Rhythmicity=r;} + double GetRhythmicity() const {return Rhythmicity;} /** Sets wind vortex, clockwise as seen from a point in front of aircraft, looking aft. Units are radians/second. */ - inline void SetWindFromClockwise(double wC) { wind_from_clockwise=wC; } - inline double GetWindFromClockwise(void) const {return wind_from_clockwise;} + void SetWindFromClockwise(double wC) { wind_from_clockwise=wC; } + double GetWindFromClockwise(void) const {return wind_from_clockwise;} - inline double GetTurbPQR(int idx) const {return vTurbPQR(idx);} + double GetTurbPQR(int idx) const {return vTurbPQR(idx);} double GetTurbMagnitude(void) const {return Magnitude;} FGColumnVector3& GetTurbDirection(void) {return vDirection;} - inline FGColumnVector3& GetTurbPQR(void) {return vTurbPQR;} + FGColumnVector3& GetTurbPQR(void) {return vTurbPQR;} protected: double rho; @@ -217,16 +251,15 @@ protected: FGColumnVector3 vDirectiondAccelDt; FGColumnVector3 vDirectionAccel; FGColumnVector3 vDirection; - FGColumnVector3 vTurbulence; FGColumnVector3 vTurbulenceGrad; FGColumnVector3 vBodyTurbGrad; FGColumnVector3 vTurbPQR; - FGColumnVector3 vWindNED; double psiw; - + FGColumnVector3 vTotalWindNED; + FGColumnVector3 vWindNED; FGColumnVector3 vGustNED; - bool bgustSet; + FGColumnVector3 vTurbulenceNED; /// Calculate the atmosphere for the given altitude, including effects of temperature deviation. void Calculate(double altitude); diff --git a/src/FDM/JSBSim/models/FGAuxiliary.cpp b/src/FDM/JSBSim/models/FGAuxiliary.cpp index fdde24a97..a04ff725d 100755 --- a/src/FDM/JSBSim/models/FGAuxiliary.cpp +++ b/src/FDM/JSBSim/models/FGAuxiliary.cpp @@ -169,10 +169,10 @@ bool FGAuxiliary::Run() } else if (GroundReactions->GetWOW() && vUVW(eU) < 30) { double factor = (vUVW(eU) - 10.0)/20.0; vAeroPQR = vPQR + factor*Atmosphere->GetTurbPQR(); - vAeroUVW = vUVW + factor*Propagate->GetTl2b()*(Atmosphere->GetWindNED()+Atmosphere->GetGustNED()); + vAeroUVW = vUVW + factor*Propagate->GetTl2b()*Atmosphere->GetTotalWindNED(); } else { vAeroPQR = vPQR + Atmosphere->GetTurbPQR(); - vAeroUVW = vUVW + Propagate->GetTl2b()*(Atmosphere->GetWindNED()+Atmosphere->GetGustNED()); + vAeroUVW = vUVW + Propagate->GetTl2b()*Atmosphere->GetTotalWindNED(); } Vt = vAeroUVW.Magnitude(); @@ -291,7 +291,7 @@ double FGAuxiliary::GetHeadWind(void) const double psiw,vw; psiw = Atmosphere->GetWindPsi(); - vw = Atmosphere->GetWindNED().Magnitude(); + vw = Atmosphere->GetTotalWindNED().Magnitude(); return vw*cos(psiw - Propagate->GetEuler(ePsi)); } @@ -303,7 +303,7 @@ double FGAuxiliary::GetCrossWind(void) const double psiw,vw; psiw = Atmosphere->GetWindPsi(); - vw = Atmosphere->GetWindNED().Magnitude(); + vw = Atmosphere->GetTotalWindNED().Magnitude(); return vw*sin(psiw - Propagate->GetEuler(ePsi)); } diff --git a/src/FDM/JSBSim/models/FGFCS.cpp b/src/FDM/JSBSim/models/FGFCS.cpp index e68827351..95abe9ebb 100644 --- a/src/FDM/JSBSim/models/FGFCS.cpp +++ b/src/FDM/JSBSim/models/FGFCS.cpp @@ -557,7 +557,7 @@ bool FGFCS::Load(Element* el, SystemType systype) if (!fname.empty()) { property_element = el->FindElement("property"); - if (property_element) cout << endl << " Declared properties" << endl << endl; + if (property_element && debug_lvl > 0) cout << endl << " Declared properties" << endl << endl; while (property_element) { double value=0.0; if ( ! property_element->GetAttributeValue("value").empty()) @@ -573,7 +573,8 @@ bool FGFCS::Load(Element* el, SystemType systype) } else { interface_properties.push_back(new double(value)); PropertyManager->Tie(interface_property_string, interface_properties.back()); - cout << " " << interface_property_string << " (initial value: " << value << ")" << endl; + if (debug_lvl > 0) + cout << " " << interface_property_string << " (initial value: " << value << ")" << endl; } @@ -600,7 +601,8 @@ bool FGFCS::Load(Element* el, SystemType systype) channel_element = document->FindElement("channel"); while (channel_element) { - cout << endl << highint << fgblue << " Channel " + if (debug_lvl > 0) + cout << endl << highint << fgblue << " Channel " << normint << channel_element->GetAttributeValue("name") << reset << endl; component_element = channel_element->GetElement(); diff --git a/src/FDM/JSBSim/models/FGMassBalance.cpp b/src/FDM/JSBSim/models/FGMassBalance.cpp index bbc5805a1..ecdd86a52 100644 --- a/src/FDM/JSBSim/models/FGMassBalance.cpp +++ b/src/FDM/JSBSim/models/FGMassBalance.cpp @@ -130,7 +130,7 @@ bool FGMassBalance::Load(Element* el) element = el->FindNextElement("pointmass"); } - Weight = EmptyWeight + Propulsion->GetTanksWeight() + GetPointMassWeight() + Weight = EmptyWeight + Propulsion->GetTanksWeight() + GetTotalPointMassWeight() + BuoyantForces->GetGasMass()*slugtolb; Mass = lbtoslug*Weight; @@ -149,7 +149,7 @@ bool FGMassBalance::Run(void) if (FGModel::Run()) return true; if (FDMExec->Holding()) return false; - Weight = EmptyWeight + Propulsion->GetTanksWeight() + GetPointMassWeight() + Weight = EmptyWeight + Propulsion->GetTanksWeight() + GetTotalPointMassWeight() + BuoyantForces->GetGasMass()*slugtolb; Mass = lbtoslug*Weight; @@ -205,8 +205,6 @@ bool FGMassBalance::Run(void) void FGMassBalance::AddPointMass(Element* el) { - char tmp[80]; - Element* loc_element = el->FindElement("location"); string pointmass_name = el->GetAttributeValue("name"); if (!loc_element) { @@ -216,18 +214,14 @@ void FGMassBalance::AddPointMass(Element* el) double w = el->FindElementValueAsNumberConvertTo("weight", "LBS"); FGColumnVector3 vXYZ = loc_element->FindElementTripletConvertTo("IN"); - PointMasses.push_back(new PointMass(w, vXYZ)); - int num = PointMasses.size()-1; - - snprintf(tmp, 80, "inertia/pointmass-weight-lbs[%u]", num); - PropertyManager->Tie( tmp, this, num, &FGMassBalance::GetPointMassWeight, - &FGMassBalance::SetPointMassWeight); + PointMasses.push_back(new PointMass(w, vXYZ)); + PointMasses.back()->bind(PropertyManager, PointMasses.size()-1); } //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -double FGMassBalance::GetPointMassWeight(void) +double FGMassBalance::GetTotalPointMassWeight(void) { double PM_total_weight = 0.0; diff --git a/src/FDM/JSBSim/models/FGMassBalance.h b/src/FDM/JSBSim/models/FGMassBalance.h index b3008b283..64f09253d 100644 --- a/src/FDM/JSBSim/models/FGMassBalance.h +++ b/src/FDM/JSBSim/models/FGMassBalance.h @@ -50,6 +50,10 @@ DEFINITIONS #define ID_MASSBALANCE "$Id$" +#if defined(WIN32) && !defined(__CYGWIN__) +#define snprintf _snprintf +#endif + /*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% FORWARD DECLARATIONSS %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/ @@ -150,23 +154,13 @@ public: inline void SetBaseCG(const FGColumnVector3& CG) {vbaseXYZcg = vXYZcg = CG;} void AddPointMass(Element* el); - double GetPointMassWeight(void); - double GetPointMassWeight(int idx) const { - if (idx < (int)PointMasses.size()) return(PointMasses[idx]->Weight); - else return 0.0; - } - - void SetPointMassWeight(int idx, double pmw) { - if (idx < (int)PointMasses.size()) { - PointMasses[idx]->Weight = pmw; - } - } + double GetTotalPointMassWeight(void); FGColumnVector3& GetPointMassMoment(void); FGMatrix33& GetJ(void) {return mJ;} FGMatrix33& GetJinv(void) {return mJinv;} void SetAircraftBaseInertias(FGMatrix33 BaseJ) {baseJ = BaseJ;} - + private: double Weight; double EmptyWeight; @@ -183,12 +177,33 @@ private: FGMatrix33& CalculatePMInertias(void); struct PointMass { + char tmp[80]; PointMass(double w, FGColumnVector3& vXYZ) { Weight = w; Location = vXYZ; } FGColumnVector3 Location; double Weight; + double GetPointMassLocation(int axis) const {return Location(axis);} + void SetPointMassLocation(int axis, double value) {Location(axis) = value;} + void SetPointMassWeight(double wt) {Weight = wt;} + double GetPointMassWeight(void) const {return Weight;} + + void bind(FGPropertyManager* PropertyManager, int num) { + snprintf(tmp, 80, "inertia/pointmass-weight-lbs[%u]", num); + PropertyManager->Tie( tmp, this, &PointMass::GetPointMassWeight, + &PointMass::SetPointMassWeight); + + snprintf(tmp, 80, "inertia/pointmass-location-X-inches[%u]", num); + PropertyManager->Tie( tmp, this, eX, &PointMass::GetPointMassLocation, + &PointMass::SetPointMassLocation); + snprintf(tmp, 80, "inertia/pointmass-location-Y-inches[%u]", num); + PropertyManager->Tie( tmp, this, eY, &PointMass::GetPointMassLocation, + &PointMass::SetPointMassLocation); + snprintf(tmp, 80, "inertia/pointmass-location-Z-inches[%u]", num); + PropertyManager->Tie( tmp, this, eZ, &PointMass::GetPointMassLocation, + &PointMass::SetPointMassLocation); + } }; vector PointMasses; diff --git a/src/FDM/JSBSim/models/FGOutput.cpp b/src/FDM/JSBSim/models/FGOutput.cpp index c84433c99..032245e85 100644 --- a/src/FDM/JSBSim/models/FGOutput.cpp +++ b/src/FDM/JSBSim/models/FGOutput.cpp @@ -371,7 +371,7 @@ void FGOutput::DelimitedOutput(string fname) outstream << Atmosphere->GetPressure() << delimeter; outstream << Atmosphere->GetTurbMagnitude() << delimeter; outstream << Atmosphere->GetTurbDirection().Dump(delimeter) << delimeter; - outstream << Atmosphere->GetWindNED().Dump(delimeter); + outstream << Atmosphere->GetTotalWindNED().Dump(delimeter); } if (SubSystems & ssMassProps) { outstream << delimeter; @@ -809,7 +809,7 @@ void FGOutput::SocketOutput(void) socket->Append(Atmosphere->GetPressure()); socket->Append(Atmosphere->GetTurbMagnitude()); socket->Append(Atmosphere->GetTurbDirection().Dump(",")); - socket->Append(Atmosphere->GetWindNED().Dump(",")); + socket->Append(Atmosphere->GetTotalWindNED().Dump(",")); } if (SubSystems & ssMassProps) { socket->Append(MassBalance->GetJ()(1,1)); diff --git a/src/FDM/JSBSim/models/FGPropulsion.cpp b/src/FDM/JSBSim/models/FGPropulsion.cpp index 1ccb7c742..2ec2e2bb3 100644 --- a/src/FDM/JSBSim/models/FGPropulsion.cpp +++ b/src/FDM/JSBSim/models/FGPropulsion.cpp @@ -462,9 +462,13 @@ FGMatrix33& FGPropulsion::CalculateTankInertias(void) tankJ = FGMatrix33(); - for (unsigned int i=0; iGetPointmassInertia( lbtoslug * Tanks[i]->GetContents(), Tanks[i]->GetXYZ() ); + tankJ(1,1) += Tanks[i]->GetIxx(); + tankJ(2,2) += Tanks[i]->GetIyy(); + tankJ(3,3) += Tanks[i]->GetIzz(); + } return tankJ; } diff --git a/src/FDM/JSBSim/models/atmosphere/FGMSIS.cpp b/src/FDM/JSBSim/models/atmosphere/FGMSIS.cpp index f5079f99e..6eeff5af8 100755 --- a/src/FDM/JSBSim/models/atmosphere/FGMSIS.cpp +++ b/src/FDM/JSBSim/models/atmosphere/FGMSIS.cpp @@ -181,7 +181,7 @@ bool MSIS::Run(void) if (turbType != ttNone) { Turbulence(); - vWindNED += vTurbulence; + vWindNED += vTurbulenceNED; } if (vWindNED(1) != 0.0) psiw = atan2( vWindNED(2), vWindNED(1) ); @@ -1643,14 +1643,14 @@ void MSIS::Debug(int from) } if (debug_lvl & 32) { // Turbulence if (first_pass && from == 2) { - cout << "vTurbulence(X), vTurbulence(Y), vTurbulence(Z), " + cout << "vTurbulenceNED(X), vTurbulenceNED(Y), vTurbulenceNED(Z), " << "vTurbulenceGrad(X), vTurbulenceGrad(Y), vTurbulenceGrad(Z), " << "vDirection(X), vDirection(Y), vDirection(Z), " << "Magnitude, " << "vTurbPQR(P), vTurbPQR(Q), vTurbPQR(R), " << endl; } if (from == 2) { - cout << vTurbulence << ", " << vTurbulenceGrad << ", " << vDirection << ", " << Magnitude << ", " << vTurbPQR << endl; + cout << vTurbulenceNED << ", " << vTurbulenceGrad << ", " << vDirection << ", " << Magnitude << ", " << vTurbPQR << endl; } } if (debug_lvl & 64) { diff --git a/src/FDM/JSBSim/models/atmosphere/FGMars.cpp b/src/FDM/JSBSim/models/atmosphere/FGMars.cpp index 544484587..d110437dd 100755 --- a/src/FDM/JSBSim/models/atmosphere/FGMars.cpp +++ b/src/FDM/JSBSim/models/atmosphere/FGMars.cpp @@ -59,83 +59,10 @@ FGMars::FGMars(FGFDMExec* fdmex) : FGAtmosphere(fdmex) Name = "FGMars"; Reng = 53.5 * 44.01; -/* - lastIndex = 0; - h = 0.0; - psiw = 0.0; - - MagnitudedAccelDt = MagnitudeAccel = Magnitude = 0.0; -// turbType = ttNone; - turbType = ttStandard; -// turbType = ttBerndt; - TurbGain = 0.0; - TurbRate = 1.0; -*/ - bind(); Debug(0); } -//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -/* -FGMars::~FGMars() -{ - Debug(1); -} -*/ -//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - -bool FGMars::InitModel(void) -{ - FGModel::InitModel(); - - Calculate(h); - SLtemperature = intTemperature; - SLpressure = intPressure; - SLdensity = intDensity; - SLsoundspeed = sqrt(SHRatio*Reng*intTemperature); - rSLtemperature = 1.0/intTemperature; - rSLpressure = 1.0/intPressure; - rSLdensity = 1.0/intDensity; - rSLsoundspeed = 1.0/SLsoundspeed; - temperature = &intTemperature; - pressure = &intPressure; - density = &intDensity; - - useExternal=false; - - return true; -} - -//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - -bool FGMars::Run(void) -{ - if (FGModel::Run()) return true; - if (FDMExec->Holding()) return false; - - //do temp, pressure, and density first - if (!useExternal) { - h = Propagate->Geth(); - Calculate(h); - } - - if (turbType != ttNone) { - Turbulence(); - vWindNED += vTurbulence; - } - - if (vWindNED(1) != 0.0) psiw = atan2( vWindNED(2), vWindNED(1) ); - - if (psiw < 0) psiw += 2*M_PI; - - soundspeed = sqrt(SHRatio*Reng*(*temperature)); - - Debug(2); - - return false; -} - //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% void FGMars::Calculate(double altitude) @@ -155,122 +82,6 @@ void FGMars::Calculate(double altitude) //cout << "Atmosphere: h=" << altitude << " rho= " << intDensity << endl; } -//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - -// square a value, but preserve the original sign - -static inline double -square_signed (double value) -{ - if (value < 0) - return value * value * -1; - else - return value * value; -} - -void FGMars::Turbulence(void) -{ - switch (turbType) { - case ttStandard: { - vDirectiondAccelDt(eX) = 1 - 2.0*(double(rand())/double(RAND_MAX)); - vDirectiondAccelDt(eY) = 1 - 2.0*(double(rand())/double(RAND_MAX)); - vDirectiondAccelDt(eZ) = 1 - 2.0*(double(rand())/double(RAND_MAX)); - - MagnitudedAccelDt = 1 - 2.0*(double(rand())/double(RAND_MAX)) - Magnitude; - // Scale the magnitude so that it moves - // away from the peaks - MagnitudedAccelDt = ((MagnitudedAccelDt - Magnitude) / - (1 + fabs(Magnitude))); - MagnitudeAccel += MagnitudedAccelDt*rate*TurbRate*State->Getdt(); - Magnitude += MagnitudeAccel*rate*State->Getdt(); - - vDirectiondAccelDt.Normalize(); - - // deemphasise non-vertical forces - vDirectiondAccelDt(eX) = square_signed(vDirectiondAccelDt(eX)); - vDirectiondAccelDt(eY) = square_signed(vDirectiondAccelDt(eY)); - - vDirectionAccel += vDirectiondAccelDt*rate*TurbRate*State->Getdt(); - vDirectionAccel.Normalize(); - vDirection += vDirectionAccel*rate*State->Getdt(); - - vDirection.Normalize(); - - // Diminish turbulence within three wingspans - // of the ground - vTurbulence = TurbGain * Magnitude * vDirection; - double HOverBMAC = Auxiliary->GetHOverBMAC(); - if (HOverBMAC < 3.0) - vTurbulence *= (HOverBMAC / 3.0) * (HOverBMAC / 3.0); - - vTurbulenceGrad = TurbGain*MagnitudeAccel * vDirection; - - vBodyTurbGrad = Propagate->GetTl2b()*vTurbulenceGrad; - vTurbPQR(eP) = vBodyTurbGrad(eY)/Aircraft->GetWingSpan(); -// if (Aircraft->GetHTailArm() != 0.0) -// vTurbPQR(eQ) = vBodyTurbGrad(eZ)/Aircraft->GetHTailArm(); -// else -// vTurbPQR(eQ) = vBodyTurbGrad(eZ)/10.0; - - if (Aircraft->GetVTailArm()) - vTurbPQR(eR) = vBodyTurbGrad(eX)/Aircraft->GetVTailArm(); - else - vTurbPQR(eR) = vBodyTurbGrad(eX)/10.0; - - // Clear the horizontal forces - // actually felt by the plane, now - // that we've used them to calculate - // moments. - vTurbulence(eX) = 0.0; - vTurbulence(eY) = 0.0; - - break; - } - case ttBerndt: { - vDirectiondAccelDt(eX) = 1 - 2.0*(double(rand())/double(RAND_MAX)); - vDirectiondAccelDt(eY) = 1 - 2.0*(double(rand())/double(RAND_MAX)); - vDirectiondAccelDt(eZ) = 1 - 2.0*(double(rand())/double(RAND_MAX)); - - - MagnitudedAccelDt = 1 - 2.0*(double(rand())/double(RAND_MAX)) - Magnitude; - MagnitudeAccel += MagnitudedAccelDt*rate*State->Getdt(); - Magnitude += MagnitudeAccel*rate*State->Getdt(); - - vDirectiondAccelDt.Normalize(); - vDirectionAccel += vDirectiondAccelDt*rate*State->Getdt(); - vDirectionAccel.Normalize(); - vDirection += vDirectionAccel*rate*State->Getdt(); - - // Diminish z-vector within two wingspans - // of the ground - double HOverBMAC = Auxiliary->GetHOverBMAC(); - if (HOverBMAC < 2.0) - vDirection(eZ) *= HOverBMAC / 2.0; - - vDirection.Normalize(); - - vTurbulence = TurbGain*Magnitude * vDirection; - vTurbulenceGrad = TurbGain*MagnitudeAccel * vDirection; - - vBodyTurbGrad = Propagate->GetTl2b()*vTurbulenceGrad; - vTurbPQR(eP) = vBodyTurbGrad(eY)/Aircraft->GetWingSpan(); - if (Aircraft->GetHTailArm() != 0.0) - vTurbPQR(eQ) = vBodyTurbGrad(eZ)/Aircraft->GetHTailArm(); - else - vTurbPQR(eQ) = vBodyTurbGrad(eZ)/10.0; - - if (Aircraft->GetVTailArm()) - vTurbPQR(eR) = vBodyTurbGrad(eX)/Aircraft->GetVTailArm(); - else - vTurbPQR(eR) = vBodyTurbGrad(eX)/10.0; - - break; - } - default: - break; - } -} - //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% // The bitmasked value choices are as follows: // unset: In this case (the default) JSBSim would only print @@ -310,13 +121,13 @@ void FGMars::Debug(int from) } if (debug_lvl & 32) { // Turbulence if (first_pass && from == 2) { - cout << "vTurbulence(X), vTurbulence(Y), vTurbulence(Z), " + cout << "vTurbulenceNED(X), vTurbulenceNED(Y), vTurbulenceNED(Z), " << "vTurbulenceGrad(X), vTurbulenceGrad(Y), vTurbulenceGrad(Z), " << "vDirection(X), vDirection(Y), vDirection(Z), " << "Magnitude, " << "vTurbPQR(P), vTurbPQR(Q), vTurbPQR(R), " << endl; } else if (from == 2) { - cout << vTurbulence << ", " << vTurbulenceGrad << ", " << vDirection << ", " << Magnitude << ", " << vTurbPQR << endl; + cout << vTurbulenceNED << ", " << vTurbulenceGrad << ", " << vDirection << ", " << Magnitude << ", " << vTurbPQR << endl; } } if (debug_lvl & 64) { diff --git a/src/FDM/JSBSim/models/atmosphere/FGMars.h b/src/FDM/JSBSim/models/atmosphere/FGMars.h index 279f771b0..4841e71ec 100755 --- a/src/FDM/JSBSim/models/atmosphere/FGMars.h +++ b/src/FDM/JSBSim/models/atmosphere/FGMars.h @@ -67,111 +67,12 @@ CLASS DECLARATION class FGMars : public FGAtmosphere { public: - /// Constructor FGMars(FGFDMExec*); - /// Destructor - //~FGMars(); - /** Runs the Martian atmosphere model; called by the Executive - @return false if no error */ - bool Run(void); - - bool InitModel(void); - - /// Returns the temperature in degrees Rankine. - inline double GetTemperature(void) const {return *temperature;} - /** Returns the density in slugs/ft^3. - This function may only be used if Run() is called first. */ - inline double GetDensity(void) const {return *density;} - /// Returns the pressure in psf. - inline double GetPressure(void) const {return *pressure;} - /// Returns the speed of sound in ft/sec. - inline double GetSoundSpeed(void) const {return soundspeed;} - - /// Returns the sea level temperature in degrees Rankine. - inline double GetTemperatureSL(void) const { return SLtemperature; } - /// Returns the sea level density in slugs/ft^3 - inline double GetDensitySL(void) const { return SLdensity; } - /// Returns the sea level pressure in psf. - inline double GetPressureSL(void) const { return SLpressure; } - /// Returns the sea level speed of sound in ft/sec. - inline double GetSoundSpeedSL(void) const { return SLsoundspeed; } - - /// Returns the ratio of at-altitude temperature over the sea level value. - inline double GetTemperatureRatio(void) const { return (*temperature)*rSLtemperature; } - /// Returns the ratio of at-altitude density over the sea level value. - inline double GetDensityRatio(void) const { return (*density)*rSLdensity; } - /// Returns the ratio of at-altitude pressure over the sea level value. - inline double GetPressureRatio(void) const { return (*pressure)*rSLpressure; } - /// Returns the ratio of at-altitude sound speed over the sea level value. - inline double GetSoundSpeedRatio(void) const { return soundspeed*rSLsoundspeed; } - - /// Tells the simulator to use an externally calculated atmosphere model. - void UseExternal(void); - /// Tells the simulator to use the internal atmosphere model. - void UseInternal(void); //this is the default - /// Gets the boolean that tells if the external atmosphere model is being used. - bool External(void) { return useExternal; } - - /// Provides the external atmosphere model with an interface to set the temperature. - inline void SetExTemperature(double t) { exTemperature=t; } - /// Provides the external atmosphere model with an interface to set the density. - inline void SetExDensity(double d) { exDensity=d; } - /// Provides the external atmosphere model with an interface to set the pressure. - inline void SetExPressure(double p) { exPressure=p; } - - /// Sets the wind components in NED frame. - inline void SetWindNED(double wN, double wE, double wD) { vWindNED(1)=wN; vWindNED(2)=wE; vWindNED(3)=wD;} - - /// Retrieves the wind components in NED frame. - inline FGColumnVector3& GetWindNED(void) { return vWindNED; } - - /** Retrieves the wind direction. The direction is defined as north=0 and - increases counterclockwise. The wind heading is returned in radians.*/ - inline double GetWindPsi(void) const { return psiw; } - - inline void SetTurbGain(double tt) {TurbGain = tt;} - inline void SetTurbRate(double tt) {TurbRate = tt;} - - inline double GetTurbPQR(int idx) const {return vTurbPQR(idx);} - inline FGColumnVector3& GetTurbPQR(void) {return vTurbPQR;} - - //void bind(void); - void unbind(void); - - -private: - double rho; - - enum tType {ttStandard, ttBerndt, ttNone} turbType; - - int lastIndex; - double h; - double htab[8]; - double SLtemperature,SLdensity,SLpressure,SLsoundspeed; - double rSLtemperature,rSLdensity,rSLpressure,rSLsoundspeed; //reciprocals - double *temperature,*density,*pressure; - double soundspeed; - bool useExternal; - double exTemperature,exDensity,exPressure; - double intTemperature, intDensity, intPressure; - - double MagnitudedAccelDt, MagnitudeAccel, Magnitude; - double TurbGain; - double TurbRate; - FGColumnVector3 vDirectiondAccelDt; - FGColumnVector3 vDirectionAccel; - FGColumnVector3 vDirection; - FGColumnVector3 vTurbulence; - FGColumnVector3 vTurbulenceGrad; - FGColumnVector3 vBodyTurbGrad; - FGColumnVector3 vTurbPQR; - - FGColumnVector3 vWindNED; - double psiw; +private: void Calculate(double altitude); - void Turbulence(void); + void Debug(int from); }; diff --git a/src/FDM/JSBSim/models/propulsion/FGElectric.h b/src/FDM/JSBSim/models/propulsion/FGElectric.h index 58e5bc2ce..ed5a93f57 100644 --- a/src/FDM/JSBSim/models/propulsion/FGElectric.h +++ b/src/FDM/JSBSim/models/propulsion/FGElectric.h @@ -82,13 +82,14 @@ public: double Calculate(void); double GetPowerAvailable(void) {return PowerAvailable;} - double CalcFuelNeed(void); double getRPM(void) {return RPM;} string GetEngineLabels(string delimeter); string GetEngineValues(string delimeter); private: + double CalcFuelNeed(void); + double BrakeHorsePower; double PowerAvailable; diff --git a/src/FDM/JSBSim/models/propulsion/FGEngine.cpp b/src/FDM/JSBSim/models/propulsion/FGEngine.cpp index cededea8b..49f67e2f4 100644 --- a/src/FDM/JSBSim/models/propulsion/FGEngine.cpp +++ b/src/FDM/JSBSim/models/propulsion/FGEngine.cpp @@ -64,7 +64,7 @@ FGEngine::FGEngine(FGFDMExec* exec, Element* engine_element, int engine_number) Type = etUnknown; X = Y = Z = 0.0; EnginePitch = EngineYaw = 0.0; - SLFuelFlowMax = SLOxiFlowMax = 0.0; + SLFuelFlowMax = 0.0; MaxThrottle = 1.0; MinThrottle = 0.0; @@ -117,8 +117,11 @@ FGEngine::FGEngine(FGFDMExec* exec, Element* engine_element, int engine_number) char property_name[80]; snprintf(property_name, 80, "propulsion/engine[%d]/set-running", EngineNumber); - PropertyManager->Tie( property_name, (FGEngine*)this, &FGEngine::GetRunning, - &FGEngine::SetRunning ); + PropertyManager->Tie( property_name, this, &FGEngine::GetRunning, &FGEngine::SetRunning ); + snprintf(property_name, 80, "propulsion/engine[%u]/thrust-lbs", EngineNumber); + PropertyManager->Tie( property_name, this, &FGEngine::GetThrust); + snprintf(property_name, 80, "propulsion/engine[%u]/fuel-flow-rate-pps", EngineNumber); + PropertyManager->Tie( property_name, this, &FGEngine::GetFuelFlowRate); Debug(0); } @@ -139,7 +142,7 @@ void FGEngine::ResetToIC(void) Throttle = 0.0; Mixture = 1.0; Starter = false; - FuelNeed = OxidizerNeed = 0.0; + FuelExpended = 0.0; Starved = Running = Cranking = false; PctPower = 0.0; TrimMode = false; @@ -153,6 +156,7 @@ void FGEngine::ResetToIC(void) // derived class' Calculate() function before any other calculations are done. // This base class method removes fuel from the fuel tanks as appropriate, // and sets the starved flag if necessary. +// This version of the fuel consumption code should never see an oxidizer tank. void FGEngine::ConsumeFuel(void) { @@ -160,22 +164,20 @@ void FGEngine::ConsumeFuel(void) if (TrimMode) return; unsigned int i; - double Fshortage, Oshortage, TanksWithFuel, TanksWithOxidizer; + double Fshortage, TanksWithFuel; FGTank* Tank; - bool haveOxTanks = false; - Fshortage = Oshortage = TanksWithFuel = TanksWithOxidizer = 0.0; + Fshortage = TanksWithFuel = 0.0; // count how many assigned tanks have fuel for (i=0; iGetTank(SourceTanks[i]); if (Tank->GetType() == FGTank::ttFUEL){ if (Tank->GetContents() > 0.0) ++TanksWithFuel; - } else if (Tank->GetType() == FGTank::ttOXIDIZER) { - haveOxTanks = true; - if (Tank->GetContents() > 0.0) ++TanksWithOxidizer; + } else { + cerr << "No oxidizer tanks should be used for this engine type." << endl; } } - if (TanksWithFuel==0 || (haveOxTanks && TanksWithOxidizer==0)) { + if (TanksWithFuel==0) { Starved = true; return; } @@ -184,12 +186,12 @@ void FGEngine::ConsumeFuel(void) Tank = Propulsion->GetTank(SourceTanks[i]); if (Tank->GetType() == FGTank::ttFUEL) { Fshortage += Tank->Drain(CalcFuelNeed()/TanksWithFuel); - } else if (Tank->GetType() == FGTank::ttOXIDIZER) { - Oshortage += Tank->Drain(CalcOxidizerNeed()/TanksWithOxidizer); + } else { + cerr << "No oxidizer tanks should be used for this engine type." << endl; } } - if (Fshortage < 0.00 || Oshortage < 0.00) Starved = true; + if (Fshortage < 0.00) Starved = true; else Starved = false; } @@ -197,16 +199,10 @@ void FGEngine::ConsumeFuel(void) double FGEngine::CalcFuelNeed(void) { - FuelNeed = SLFuelFlowMax*PctPower*State->Getdt()*Propulsion->GetRate(); - return FuelNeed; -} - -//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - -double FGEngine::CalcOxidizerNeed(void) -{ - OxidizerNeed = SLOxiFlowMax*PctPower*State->Getdt()*Propulsion->GetRate(); - return OxidizerNeed; + double dT = State->Getdt()*Propulsion->GetRate(); + FuelFlowRate = SLFuelFlowMax*PctPower; + FuelExpended = FuelFlowRate*dT; + return FuelExpended; } //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% diff --git a/src/FDM/JSBSim/models/propulsion/FGEngine.h b/src/FDM/JSBSim/models/propulsion/FGEngine.h index 225c39a0f..dd8b1302b 100644 --- a/src/FDM/JSBSim/models/propulsion/FGEngine.h +++ b/src/FDM/JSBSim/models/propulsion/FGEngine.h @@ -149,7 +149,8 @@ public: virtual double getFuelFlow_gph () const {return FuelFlow_gph;} virtual double getFuelFlow_pph () const {return FuelFlow_pph;} - virtual double GetThrust(void) { return Thrust; } + virtual double GetFuelFlowRate(void) const {return FuelFlowRate;} + virtual double GetThrust(void) const { return Thrust; } virtual bool GetStarved(void) { return Starved; } virtual bool GetRunning(void) const { return Running; } virtual bool GetCranking(void) { return Cranking; } @@ -173,25 +174,6 @@ public: @return Thrust in pounds */ virtual double Calculate(void) {return 0.0;} - /** Reduces the fuel in the active tanks by the amount required. - This function should be called from within the - derived class' Calculate() function before any other calculations are - done. This base class method removes fuel from the fuel tanks as - appropriate, and sets the starved flag if necessary. */ - virtual void ConsumeFuel(void); - - /** The fuel need is calculated based on power levels and flow rate for that - power level. It is also turned from a rate into an actual amount (pounds) - by multiplying it by the delta T and the rate. - @return Total fuel requirement for this engine in pounds. */ - virtual double CalcFuelNeed(void); - - /** The oxidizer need is calculated based on power levels and flow rate for that - power level. It is also turned from a rate into an actual amount (pounds) - by multiplying it by the delta T and the rate. - @return Total oxidizer requirement for this engine in pounds. */ - virtual double CalcOxidizerNeed(void); - /// Sets engine placement information virtual void SetPlacement(FGColumnVector3& location, FGColumnVector3& orientation); @@ -210,6 +192,19 @@ public: virtual string GetEngineValues(string delimeter) = 0; protected: + /** Reduces the fuel in the active tanks by the amount required. + This function should be called from within the + derived class' Calculate() function before any other calculations are + done. This base class method removes fuel from the fuel tanks as + appropriate, and sets the starved flag if necessary. */ + virtual void ConsumeFuel(void); + + /** The fuel need is calculated based on power levels and flow rate for that + power level. It is also turned from a rate into an actual amount (pounds) + by multiplying it by the delta T and the rate. + @return Total fuel requirement for this engine in pounds. */ + virtual double CalcFuelNeed(void); + FGPropertyManager* PropertyManager; string Name; const int EngineNumber; @@ -218,15 +213,14 @@ protected: double EnginePitch; double EngineYaw; double SLFuelFlowMax; - double SLOxiFlowMax; double MaxThrottle; double MinThrottle; double Thrust; double Throttle; double Mixture; - double FuelNeed; - double OxidizerNeed; + double FuelExpended; + double FuelFlowRate; double PctPower; bool Starter; bool Starved; diff --git a/src/FDM/JSBSim/models/propulsion/FGNozzle.cpp b/src/FDM/JSBSim/models/propulsion/FGNozzle.cpp index 476f4a988..dc397cd7f 100644 --- a/src/FDM/JSBSim/models/propulsion/FGNozzle.cpp +++ b/src/FDM/JSBSim/models/propulsion/FGNozzle.cpp @@ -53,6 +53,12 @@ CLASS IMPLEMENTATION FGNozzle::FGNozzle(FGFDMExec* FDMExec, Element* nozzle_element, int num) : FGThruster(FDMExec, nozzle_element, num) { + if (nozzle_element->FindElement("area")) + Area = nozzle_element->FindElementValueAsNumberConvertTo("area", "FT2"); + else { + cerr << "Fatal Error: Nozzle exit area must be given in nozzle config file." << endl; + exit(-1); + } if (nozzle_element->FindElement("pe")) PE = nozzle_element->FindElementValueAsNumberConvertTo("pe", "PSF"); @@ -60,29 +66,9 @@ FGNozzle::FGNozzle(FGFDMExec* FDMExec, Element* nozzle_element, int num) cerr << "Fatal Error: Nozzle exit pressure must be given in nozzle config file." << endl; exit(-1); } - if (nozzle_element->FindElement("expr")) - ExpR = nozzle_element->FindElementValueAsNumber("expr"); - else { - cerr << "Fatal Error: Nozzle expansion ratio must be given in nozzle config file." << endl; - exit(-1); - } - if (nozzle_element->FindElement("nzl_eff")) - nzlEff = nozzle_element->FindElementValueAsNumber("nzl_eff"); - else { - cerr << "Fatal Error: Nozzle efficiency must be given in nozzle config file." << endl; - exit(-1); - } - if (nozzle_element->FindElement("diam")) - Diameter = nozzle_element->FindElementValueAsNumberConvertTo("diam", "FT"); - else { - cerr << "Fatal Error: Nozzle diameter must be given in nozzle config file." << endl; - exit(-1); - } Thrust = 0; Type = ttNozzle; - Area2 = (Diameter*Diameter/4.0)*M_PI; - AreaT = Area2/ExpR; Debug(0); } @@ -96,30 +82,18 @@ FGNozzle::~FGNozzle() //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -double FGNozzle::Calculate(double CfPc) +double FGNozzle::Calculate(double vacThrust) { double pAtm = fdmex->GetAtmosphere()->GetPressure(); - if (CfPc > 0) - Thrust = max((double)0.0, (CfPc * AreaT + (PE - pAtm)*Area2) * nzlEff); - else - Thrust = 0.0; + Thrust = max((double)0.0, vacThrust - pAtm*Area); vFn(1) = Thrust * cos(ReverserAngle); - ThrustCoeff = max((double)0.0, CfPc / ((pAtm - PE) * Area2)); - return Thrust; } //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -double FGNozzle::GetPowerRequired(void) -{ - return PE; -} - -//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - string FGNozzle::GetThrusterLabels(int id, string delimeter) { std::ostringstream buf; @@ -166,10 +140,7 @@ void FGNozzle::Debug(int from) if (debug_lvl & 1) { // Standard console startup message output if (from == 0) { // Constructor cout << " Nozzle Name: " << Name << endl; - cout << " Nozzle Exit Pressure = " << PE << endl; - cout << " Nozzle Expansion Ratio = " << ExpR << endl; - cout << " Nozzle Efficiency = " << nzlEff << endl; - cout << " Nozzle Diameter = " << Diameter << endl; + cout << " Nozzle Exit Area = " << Area << endl; } } if (debug_lvl & 2 ) { // Instantiation/Destruction notification diff --git a/src/FDM/JSBSim/models/propulsion/FGNozzle.h b/src/FDM/JSBSim/models/propulsion/FGNozzle.h index dfa3ffd7d..c9549126b 100644 --- a/src/FDM/JSBSim/models/propulsion/FGNozzle.h +++ b/src/FDM/JSBSim/models/propulsion/FGNozzle.h @@ -63,18 +63,14 @@ CLASS DOCUMENTATION @code {number} - {number} - {number} - {number} + {number} @endcode

Configuration parameters are:

-    pe -      Nozzle exit pressure.
-    expr -    Nozzle expansion ratio, Ae/At, sqft. dimensionless ratio.
-    nzl_eff - Nozzle efficiency, 0.0 - 1.0.
-    diam -    Nozzle diameter.
+    pe -      Nozzle design exit pressure.
+    area -    Nozzle area at the exit plane.
All parameters MUST be specified. @@ -94,18 +90,13 @@ public: /// Destructor ~FGNozzle(); - double Calculate(double CfPc); - double GetPowerRequired(void); + double Calculate(double vacThrust); string GetThrusterLabels(int id, string delimeter); string GetThrusterValues(int id, string delimeter); private: double PE; - double ExpR; - double nzlEff; - double Diameter; - double AreaT; - double Area2; + double Area; void Debug(int from); }; } diff --git a/src/FDM/JSBSim/models/propulsion/FGPiston.cpp b/src/FDM/JSBSim/models/propulsion/FGPiston.cpp index 7f981ce33..e23d981a8 100644 --- a/src/FDM/JSBSim/models/propulsion/FGPiston.cpp +++ b/src/FDM/JSBSim/models/propulsion/FGPiston.cpp @@ -82,6 +82,9 @@ FGPiston::FGPiston(FGFDMExec* exec, Element* el, int engine_number) MaxManifoldPressure_inHg = 28.5; BSFC = -1; + // Initialisation + volumetric_efficiency = 0.8; // Actually f(speed, load) but this will get us running + // These are internal program variables crank_counter = 0; @@ -113,8 +116,6 @@ FGPiston::FGPiston(FGFDMExec* exec, Element* el, int engine_number) BoostSwitchAltitude[i] = 0.0; BoostSwitchPressure[i] = 0.0; } - // Initialisation - volumetric_efficiency = 0.8; // Actually f(speed, load) but this will get us running // First column is thi, second is neta (combustion efficiency) Lookup_Combustion_Efficiency = new FGTable(12); @@ -202,7 +203,9 @@ Manifold_Pressure_Lookup = new if (el->FindElement("minthrottle")) MinThrottle = el->FindElementValueAsNumber("minthrottle"); if (el->FindElement("bsfc")) - BSFC = el->FindElementValueAsNumber("bsfc"); + BSFC = el->FindElementValueAsNumberConvertTo("bsfc", "LBS/HP*HR"); + if (el->FindElement("volumetric-efficiency")) + volumetric_efficiency = el->FindElementValueAsNumber("volumetric-efficiency"); if (el->FindElement("numboostspeeds")) { // Turbo- and super-charging parameters BoostSpeeds = (int)el->FindElementValueAsNumber("numboostspeeds"); if (el->FindElement("boostoverride")) @@ -236,16 +239,16 @@ Manifold_Pressure_Lookup = new } // Create a BSFC to match the engine if not provided - // The 0.8 in the equation below is volumetric efficiency if (BSFC < 0) { - BSFC = ( Displacement * MaxRPM * 0.8 ) / (9411 * MaxHP); + BSFC = ( Displacement * MaxRPM * volumetric_efficiency ) / (9411 * MaxHP); } char property_name[80]; - snprintf(property_name, 80, "propulsion/engine[%d]/power_hp", EngineNumber); + snprintf(property_name, 80, "propulsion/engine[%d]/power-hp", EngineNumber); PropertyManager->Tie(property_name, &HP); - snprintf(property_name, 80, "propulsion/engine[%d]/bsfc", EngineNumber); + snprintf(property_name, 80, "propulsion/engine[%d]/bsfc-lbs_hphr", EngineNumber); PropertyManager->Tie(property_name, &BSFC); - + snprintf(property_name, 80, "propulsion/engine[%d]/volumetric-efficiency", EngineNumber); + PropertyManager->Tie(property_name, &volumetric_efficiency); minMAP = MinManifoldPressure_inHg * inhgtopa; // inHg to Pa maxMAP = MaxManifoldPressure_inHg * inhgtopa; StarterHP = sqrt(MaxHP) * 0.4; @@ -315,7 +318,7 @@ FGPiston::~FGPiston() void FGPiston::ResetToIC(void) { FGEngine::ResetToIC(); - + ManifoldPressure_inHg = Atmosphere->GetPressure() * psftoinhg; // psf to in Hg MAP = Atmosphere->GetPressure() * psftopa; double airTemperature_degK = RankineToKelvin(Atmosphere->GetTemperature()); @@ -325,6 +328,7 @@ void FGPiston::ResetToIC(void) EGT_degC = ExhaustGasTemp_degK - 273; Thruster->SetRPM(0.0); RPM = 0.0; + OilPressure_psi = 0.0; } //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% @@ -385,7 +389,11 @@ if(HP<0.1250) double FGPiston::CalcFuelNeed(void) { - return FuelFlow_gph / 3600 * 6 * State->Getdt() * Propulsion->GetRate(); + double dT = State->Getdt() * Propulsion->GetRate(); + FuelFlow_pph = FuelFlow_gph * 6.0; // Assumes 6 lbs / gallon + FuelFlowRate = FuelFlow_pph / 3600.0; + FuelExpended = FuelFlowRate * dT; + return FuelExpended; } //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% @@ -541,7 +549,9 @@ void FGPiston::doMAP(void) } } // Boost the manifold pressure. - MAP += MAP * BoostMul[BoostSpeed] * suction_loss * RPM/RatedRPM[BoostSpeed]; + double boost_factor = BoostMul[BoostSpeed] * suction_loss * RPM/RatedRPM[BoostSpeed]; + if (boost_factor < 1.0) boost_factor = 1.0; // boost will never reduce the MAP + MAP *= boost_factor; // Now clip the manifold pressure to BCV or Wastegate setting. if(bTakeoffPos) { if(MAP > TakeoffMAP[BoostSpeed]) { @@ -574,8 +584,7 @@ void FGPiston::doMAP(void) void FGPiston::doAirFlow(void) { - -rho_air = p_amb / (R_air * T_amb); + rho_air = p_amb / (R_air * T_amb); double displacement_SI = Displacement * in3tom3; double swept_volume = (displacement_SI * (RPM/60)) / 2; double v_dot_air = swept_volume * volumetric_efficiency * suction_loss; @@ -596,9 +605,10 @@ rho_air = p_amb / (R_air * T_amb); void FGPiston::doFuelFlow(void) { double thi_sea_level = 1.3 * Mixture; // Allows an AFR of infinity:1 to 11.3075:1 - equivalence_ratio = thi_sea_level; // * p_amb_sea_level / p_amb; - double AFR = 10+(12*(1-Mixture));// mixture 10:1 to 22:1 - m_dot_fuel = m_dot_air / AFR; + equivalence_ratio = thi_sea_level * 101325.0 / p_amb; +// double AFR = 10+(12*(1-Mixture));// mixture 10:1 to 22:1 +// m_dot_fuel = m_dot_air / AFR; + m_dot_fuel = (m_dot_air * equivalence_ratio) / 14.7; FuelFlow_gph = m_dot_fuel * 3600 // seconds to hours * 2.2046 // kg to lb @@ -703,7 +713,7 @@ void FGPiston::doEGT(void) * Calculate the cylinder head temperature. * * Inputs: T_amb, IAS, rho_air, m_dot_fuel, calorific_value_fuel, - * combustion_efficiency, RPM + * combustion_efficiency, RPM, MaxRPM * * Outputs: CylinderHeadTemp_degK */ @@ -712,7 +722,7 @@ void FGPiston::doCHT(void) { double h1 = -95.0; double h2 = -3.95; - double h3 = -0.05; + double h3 = -140.0; // -0.05 * 2800 (default maxrpm) double arbitary_area = 1.0; double CpCylinderHead = 800.0; @@ -725,7 +735,7 @@ void FGPiston::doCHT(void) double dqdt_from_combustion = m_dot_fuel * calorific_value_fuel * combustion_efficiency * 0.33; double dqdt_forced = (h2 * m_dot_cooling_air * temperature_difference) + - (h3 * RPM * temperature_difference); + (h3 * RPM * temperature_difference / MaxRPM); double dqdt_free = h1 * temperature_difference; double dqdt_cylinder_head = dqdt_from_combustion + dqdt_forced + dqdt_free; @@ -739,7 +749,7 @@ void FGPiston::doCHT(void) /** * Calculate the oil temperature. * - * Inputs: Percentage_Power, running flag. + * Inputs: CylinderHeadTemp_degK, T_amb, OilPressure_psi. * * Outputs: OilTemp_degK */ @@ -749,15 +759,18 @@ void FGPiston::doOilTemperature(void) double idle_percentage_power = 0.023; // approximately double target_oil_temp; // Steady state oil temp at the current engine conditions double time_constant; // The time constant for the differential equation + double efficiency = 0.667; // The aproximate oil cooling system efficiency // FIXME: may vary by engine - 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 - } +// Target oil temp is interpolated between ambient temperature and Cylinder Head Tempurature +// target_oil_temp = ( T_amb * efficiency ) + (CylinderHeadTemp_degK *(1-efficiency)) ; + target_oil_temp = CylinderHeadTemp_degK + efficiency * (T_amb - CylinderHeadTemp_degK) ; + + if (OilPressure_psi > 5.0 ) { + time_constant = 5000 / OilPressure_psi; // Guess at a time constant for circulated oil. + // The higher the pressure the faster it reaches + // target temperature. Oil pressure should be about + // 60 PSI yielding a TC of about 80. } else { - target_oil_temp = RankineToKelvin(Atmosphere->GetTemperature()); time_constant = 1000; // Time constant for engine-off; reflects the fact // that oil is no longer getting circulated } @@ -771,7 +784,7 @@ void FGPiston::doOilTemperature(void) /** * Calculate the oil pressure. * - * Inputs: RPM + * Inputs: RPM, MaxRPM, OilTemp_degK * * Outputs: OilPressure_psi */ diff --git a/src/FDM/JSBSim/models/propulsion/FGPiston.h b/src/FDM/JSBSim/models/propulsion/FGPiston.h index 14d86e3b0..d4fa080f3 100644 --- a/src/FDM/JSBSim/models/propulsion/FGPiston.h +++ b/src/FDM/JSBSim/models/propulsion/FGPiston.h @@ -27,7 +27,7 @@ HISTORY -------------------------------------------------------------------------------- 09/12/2000 JSB Created 10/01/2001 DPM Modified to use equations from Dave Luff's piston model. - +11/01/2008 RKJ Modified piston engine model for more general use. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% SENTRY %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/ @@ -66,15 +66,19 @@ CLASS DOCUMENTATION @code - {number} - {number} + {number} + {number} {number} + {number} {number} {number} {number} + {number} {number} {number} {number} + {number} + {number} {0 | 1} {number} {number} @@ -164,6 +168,7 @@ CLASS DOCUMENTATION @author Jon S. Berndt (Engine framework code and framework-related mods) @author Dave Luff (engine operational code) @author David Megginson (initial porting and additional code) + @author Ron Jensen (additional engine code) @version $Id$ */ @@ -279,7 +284,7 @@ private: double minMAP; // Pa double maxMAP; // Pa double MAP; // Pa - double BSFC; // unitless + double BSFC; // brake specific fuel consumption [lbs/horsepower*hour // // Inputs (in addition to those in FGEngine). diff --git a/src/FDM/JSBSim/models/propulsion/FGRocket.cpp b/src/FDM/JSBSim/models/propulsion/FGRocket.cpp index ada946f9d..66390010c 100644 --- a/src/FDM/JSBSim/models/propulsion/FGRocket.cpp +++ b/src/FDM/JSBSim/models/propulsion/FGRocket.cpp @@ -57,18 +57,20 @@ FGRocket::FGRocket(FGFDMExec* exec, Element *el, int engine_number) Element* thrust_table_element = 0; ThrustTable = 0L; BurnTime = 0.0; + previousFuelNeedPerTank = 0.0; + previousOxiNeedPerTank = 0.0; + PropellantFlowRate = 0.0; + FuelFlowRate = 0.0; + OxidizerFlowRate = 0.0; + SLOxiFlowMax = 0.0; + It = 0.0; // Defaults - Variance = 0.0; MinThrottle = 0.0; MaxThrottle = 1.0; - if (el->FindElement("shr")) - SHR = el->FindElementValueAsNumber("shr"); - if (el->FindElement("max_pc")) - maxPC = el->FindElementValueAsNumberConvertTo("max_pc", "PSF"); - if (el->FindElement("prop_eff")) - propEff = el->FindElementValueAsNumber("prop_eff"); + if (el->FindElement("isp")) + Isp = el->FindElementValueAsNumber("isp"); if (el->FindElement("maxthrottle")) MaxThrottle = el->FindElementValueAsNumber("maxthrottle"); if (el->FindElement("minthrottle")) @@ -77,21 +79,19 @@ FGRocket::FGRocket(FGFDMExec* exec, Element *el, int engine_number) SLFuelFlowMax = el->FindElementValueAsNumberConvertTo("slfuelflowmax", "LBS/SEC"); if (el->FindElement("sloxiflowmax")) SLOxiFlowMax = el->FindElementValueAsNumberConvertTo("sloxiflowmax", "LBS/SEC"); - if (el->FindElement("variance")) - Variance = el->FindElementValueAsNumber("variance"); thrust_table_element = el->FindElement("thrust_table"); if (thrust_table_element) { ThrustTable = new FGTable(PropertyManager, thrust_table_element); } + bindmodel(); + Debug(0); Type = etRocket; Flameout = false; - PC = 0.0; - kFactor = (2.0*SHR*SHR/(SHR-1.0))*pow(2.0/(SHR+1), (SHR+1)/(SHR-1)); } //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% @@ -105,41 +105,143 @@ FGRocket::~FGRocket(void) double FGRocket::Calculate(void) { - double Cf=0; + double dT = State->Getdt()*Propulsion->GetRate(); if (!Flameout && !Starved) ConsumeFuel(); + PropellantFlowRate = (FuelExpended + OxidizerExpended)/dT; Throttle = FCS->GetThrottlePos(EngineNumber); - // If there is a thrust table, it is a function of elapsed burn time. The engine - // is started when the throttle is advanced to 1.0. After that, it burns - // without regard to throttle setting. The table returns a value between zero - // and one, representing the percentage of maximum vacuum thrust being applied. + // If there is a thrust table, it is a function of propellant remaining. The + // engine is started when the throttle is advanced to 1.0. After that, it + // burns without regard to throttle setting. The table returns a value between + // zero and one, representing the percentage of maximum vacuum thrust being + // applied. + + if (ThrustTable != 0L) { // Thrust table given -> Solid fuel used - if (ThrustTable != 0L) { - if (Throttle == 1 || BurnTime > 0.0) { + if ((Throttle == 1 || BurnTime > 0.0 ) && !Starved) { BurnTime += State->Getdt(); + double TotalEngineFuelAvailable=0.0; + for (int i=0; i<(int)SourceTanks.size(); i++) + TotalEngineFuelAvailable += Propulsion->GetTank(SourceTanks[i])->GetContents(); + + VacThrust = ThrustTable->GetValue(TotalEngineFuelAvailable); + } else { + VacThrust = 0.0; + } + + } else { // liquid fueled rocket assumed + + if (Throttle < MinThrottle || Starved) { // Combustion not supported + + PctPower = Thrust = 0.0; // desired thrust + Flameout = true; + VacThrust = 0.0; + + } else { // Calculate thrust + + PctPower = Throttle / MaxThrottle; // Min and MaxThrottle range from 0.0 to 1.0, normally. + Flameout = false; + VacThrust = Isp * PropellantFlowRate; + + } + + } // End thrust calculations + + Thrust = Thruster->Calculate(VacThrust); + It += Thrust * dT; + + return Thrust; +} + +//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +// This overrides the base class ConsumeFuel() function, for special rocket +// engine processing. + +void FGRocket::ConsumeFuel(void) +{ + unsigned int i; + FGTank* Tank; + bool haveOxTanks = false; + double Fshortage=0, Oshortage=0, TanksWithFuel=0, TanksWithOxidizer=0; + + if (FuelFreeze) return; + if (TrimMode) return; + + // Count how many assigned tanks have fuel for this engine at this time. + // If there is/are fuel tanks but no oxidizer tanks, this indicates + // a solid rocket is being modeled. + + for (i=0; iGetTank(SourceTanks[i]); + switch(Tank->GetType()) { + case FGTank::ttFUEL: + if (Tank->GetContents() > 0.0 && Tank->GetSelected()) ++TanksWithFuel; + break; + case FGTank::ttOXIDIZER: + haveOxTanks = true; + if (Tank->GetContents() > 0.0 && Tank->GetSelected()) ++TanksWithOxidizer; + break; + } + } + + // If this engine has burned out, it is starved. + + if (TanksWithFuel==0 || (haveOxTanks && TanksWithOxidizer==0)) { + Starved = true; + return; + } + + // Expend fuel from the engine's tanks if the tank is selected as a source + // for this engine. + + double fuelNeedPerTank = CalcFuelNeed()/TanksWithFuel; + double oxiNeedPerTank = CalcOxidizerNeed()/TanksWithOxidizer; + + for (i=0; iGetTank(SourceTanks[i]); + if ( ! Tank->GetSelected()) continue; // If this tank is not selected as a source, skip it. + switch(Tank->GetType()) { + case FGTank::ttFUEL: + Fshortage += Tank->Drain(2.0*fuelNeedPerTank - previousFuelNeedPerTank); + previousFuelNeedPerTank = fuelNeedPerTank; + break; + case FGTank::ttOXIDIZER: + Oshortage += Tank->Drain(2.0*oxiNeedPerTank - previousOxiNeedPerTank); + previousOxiNeedPerTank = oxiNeedPerTank; + break; } - Throttle = ThrustTable->GetValue(BurnTime); } - if (Throttle < MinThrottle || Starved) { - PctPower = Thrust = 0.0; // desired thrust - Flameout = true; - PC = 0.0; + if (Fshortage < 0.00 || (haveOxTanks && Oshortage < 0.00)) Starved = true; + else Starved = false; +} + +//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +double FGRocket::CalcFuelNeed(void) +{ + double dT = State->Getdt()*Propulsion->GetRate(); + + if (ThrustTable != 0L) { // Thrust table given - infers solid fuel + FuelFlowRate = VacThrust/Isp; // This calculates wdot (weight flow rate in lbs/sec) } else { - PctPower = Throttle / MaxThrottle; - //todo: remove Variance? - PC = maxPC*PctPower * (1.0 + Variance * ((double)rand()/(double)RAND_MAX - 0.5)); - // The Cf (below) is CF from Eqn. 3-30, "Rocket Propulsion Elements", Fifth Edition, - // George P. Sutton. Note that the thruster function GetPowerRequired() might - // be better called GetResistance() or something; this function returns the - // nozzle exit pressure. - Cf = sqrt(kFactor*(1 - pow(Thruster->GetPowerRequired()/(PC), (SHR-1)/SHR))); - Flameout = false; + FuelFlowRate = SLFuelFlowMax*PctPower; } - return Thruster->Calculate(Cf*maxPC*PctPower*propEff); + FuelExpended = FuelFlowRate*dT; // For this time step ... + return FuelExpended; +} + +//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +double FGRocket::CalcOxidizerNeed(void) +{ + double dT = State->Getdt()*Propulsion->GetRate(); + OxidizerFlowRate = SLOxiFlowMax*PctPower; + OxidizerExpended = OxidizerFlowRate*dT; + return OxidizerExpended; } //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% @@ -148,7 +250,7 @@ string FGRocket::GetEngineLabels(string delimeter) { std::ostringstream buf; - buf << Name << " Chamber Pressure (engine " << EngineNumber << " in psf)" << delimeter + buf << Name << " Total Impulse (engine " << EngineNumber << " in psf)" << delimeter << Thruster->GetThrusterLabels(EngineNumber, delimeter); return buf.str(); @@ -160,11 +262,27 @@ string FGRocket::GetEngineValues(string delimeter) { std::ostringstream buf; - buf << PC << delimeter << Thruster->GetThrusterValues(EngineNumber, delimeter); + buf << It << delimeter << Thruster->GetThrusterValues(EngineNumber, delimeter); return buf.str(); } +//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +// This funciton should tie properties to rocket engine specific properties +// that are not bound in the base class (FGEngine) code. +// +void FGRocket::bindmodel() +{ + char property_name[80]; + + snprintf(property_name, 80, "propulsion/engine[%u]/total-impulse", EngineNumber); + PropertyManager->Tie( property_name, this, &FGRocket::GetTotalImpulse); + snprintf(property_name, 80, "propulsion/engine[%u]/oxi-flow-rate-pps", EngineNumber); + PropertyManager->Tie( property_name, this, &FGRocket::GetOxiFlowRate); + snprintf(property_name, 80, "propulsion/engine[%u]/vacuum-thrust_lbs", EngineNumber); + PropertyManager->Tie( property_name, this, &FGRocket::GetVacThrust); +} + //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% // The bitmasked value choices are as follows: // unset: In this case (the default) JSBSim would only print @@ -191,14 +309,12 @@ void FGRocket::Debug(int from) if (debug_lvl & 1) { // Standard console startup message output if (from == 0) { // Constructor cout << " Engine Name: " << Name << endl; - cout << " Specific Heat Ratio = " << SHR << endl; - cout << " Maximum Chamber Pressure = " << maxPC << endl; - cout << " Propulsive Efficiency = " << propEff << endl; - cout << " MaxThrottle = " << MaxThrottle << endl; - cout << " MinThrottle = " << MinThrottle << endl; - cout << " FuelFlowMax = " << SLFuelFlowMax << endl; - cout << " OxiFlowMax = " << SLOxiFlowMax << endl; - cout << " Variance = " << Variance << endl; + cout << " Vacuum Isp = " << Isp << endl; + cout << " Maximum Throttle = " << MaxThrottle << endl; + cout << " Minimum Throttle = " << MinThrottle << endl; + cout << " Fuel Flow (max) = " << SLFuelFlowMax << endl; + cout << " Oxidizer Flow (max) = " << SLOxiFlowMax << endl; + cout << " Mixture ratio = " << SLOxiFlowMax/SLFuelFlowMax << endl; } } if (debug_lvl & 2 ) { // Instantiation/Destruction notification diff --git a/src/FDM/JSBSim/models/propulsion/FGRocket.h b/src/FDM/JSBSim/models/propulsion/FGRocket.h index 3e4479406..afefb5884 100644 --- a/src/FDM/JSBSim/models/propulsion/FGRocket.h +++ b/src/FDM/JSBSim/models/propulsion/FGRocket.h @@ -61,10 +61,7 @@ CLASS DOCUMENTATION /** Models a generic rocket engine. The rocket engine is modeled given the following parameters:
    -
  • Chamber pressure (in psf)
    • -
  • Specific heat ratio (usually about 1.2 for hydrocarbon fuel and LOX)
    • -
  • Propulsive efficiency (in percent, from 0 to 1.0)
    • -
  • Variance (in percent, from 0 to 1.0, nominally 0.05)
    • +
  • Specific Impulse (in sec)
Additionally, the following control inputs, operating characteristics, and location are required, as with all other engine types: @@ -78,12 +75,8 @@ CLASS DOCUMENTATION
  • Pitch and Yaw
    • The nozzle exit pressure (p2) is returned via a - call to FGNozzle::GetPowerRequired(). This exit pressure is used, - along with chamber pressure and specific heat ratio, to get the - thrust coefficient for the throttle setting. This thrust - coefficient is multiplied by the chamber pressure and then passed - to the nozzle Calculate() routine, where the thrust force is - determined. + call to FGNozzle::GetPowerRequired(). This exit pressure is used + to get the at-altitude thrust level. One can model the thrust of a solid rocket by providing a normalized thrust table as a function of time. For instance, the space shuttle solid rocket booster @@ -150,30 +143,62 @@ public: /** Destructor */ ~FGRocket(void); - /** Determines the thrust coefficient. - @return thrust coefficient times chamber pressure */ + /** Determines the thrust. + @return thrust */ double Calculate(void); - /** Gets the chamber pressure. - @return chamber pressure in psf. */ - double GetChamberPressure(void) {return PC;} + /** Gets the total impulse of the rocket. + @return The cumulative total impulse of the rocket up to this time.*/ + double GetTotalImpulse(void) const {return It;} /** Gets the flame-out status. The engine will "flame out" if the throttle is set below the minimum - sustainable setting. + sustainable-thrust setting. @return true if engine has flamed out. */ bool GetFlameout(void) {return Flameout;} + + double GetOxiFlowRate(void) const {return OxidizerFlowRate;} + string GetEngineLabels(string delimeter); string GetEngineValues(string delimeter); private: - double SHR; - double maxPC; - double propEff; - double kFactor; - double Variance; - double PC; + /** Reduces the fuel in the active tanks by the amount required. + This function should be called from within the + derived class' Calculate() function before any other calculations are + done. This base class method removes fuel from the fuel tanks as + appropriate, and sets the starved flag if necessary. */ + void ConsumeFuel(void); + + /** The fuel need is calculated based on power levels and flow rate for that + power level. It is also turned from a rate into an actual amount (pounds) + by multiplying it by the delta T and the rate. + @return Total fuel requirement for this engine in pounds. */ + double CalcFuelNeed(void); + + /** The oxidizer need is calculated based on power levels and flow rate for that + power level. It is also turned from a rate into an actual amount (pounds) + by multiplying it by the delta T and the rate. + @return Total oxidizer requirement for this engine in pounds. */ + double CalcOxidizerNeed(void); + + /** Returns the vacuum thrust. + @return The vacuum thrust in lbs. */ + double GetVacThrust(void) const {return VacThrust;} + + void bindmodel(void); + + double Isp; // Vacuum Isp + double It; + double MxR; // Mixture Ratio double BurnTime; + double VacThrust; + double previousFuelNeedPerTank; + double previousOxiNeedPerTank; + double OxidizerExpended; + double SLOxiFlowMax; + double OxidizerFlowRate; + double PropellantFlowRate; bool Flameout; FGTable* ThrustTable; diff --git a/src/FDM/JSBSim/models/propulsion/FGTank.cpp b/src/FDM/JSBSim/models/propulsion/FGTank.cpp index 9827c5606..5f99f402e 100644 --- a/src/FDM/JSBSim/models/propulsion/FGTank.cpp +++ b/src/FDM/JSBSim/models/propulsion/FGTank.cpp @@ -56,10 +56,11 @@ FGTank::FGTank(FGFDMExec* exec, Element* el, int tank_number) { string token; Element* element; + Element* element_Grain; Area = 1.0; Temperature = -9999.0; Auxiliary = exec->GetAuxiliary(); - Radius = Capacity = Contents = Standpipe = 0.0; + Radius = Capacity = Contents = Standpipe = Length = InnerRadius = 0.0; PropertyManager = exec->GetPropertyManager(); vXYZ.InitMatrix(); vXYZ_drain.InitMatrix(); @@ -100,6 +101,44 @@ FGTank::FGTank(FGFDMExec* exec, Element* el, int tank_number) PctFull = 0; } + // Check whether this is a solid propellant "tank". Initialize it if true. + + grainType = gtUNKNOWN; // This is the default + + element_Grain = el->FindElement("grain_config"); + if (element_Grain) { + + strGType = element_Grain->GetAttributeValue("type"); + if (strGType == "CYLINDRICAL") grainType = gtCYLINDRICAL; + else if (strGType == "ENDBURNING") grainType = gtENDBURNING; + else cerr << "Unknown propellant grain type specified" << endl; + + if (element_Grain->FindElement("length")) + Length = element_Grain->FindElementValueAsNumberConvertTo("length", "IN"); + if (element_Grain->FindElement("bore_diameter")) + InnerRadius = element_Grain->FindElementValueAsNumberConvertTo("bore_diameter", "IN")/2.0; + + // Initialize solid propellant values for debug and runtime use. + + switch (grainType) { + case gtCYLINDRICAL: + if (Radius <= InnerRadius) { + cerr << "The bore diameter should be smaller than the total grain diameter!" << endl; + exit(-1); + } + Volume = M_PI * Length * (Radius*Radius - InnerRadius*InnerRadius); // cubic inches + break; + case gtENDBURNING: + Volume = M_PI * Length * Radius * Radius; // cubic inches + break; + case gtUNKNOWN: + cerr << "Unknown grain type found in this rocket engine definition." << endl; + exit(-1); + } + Density = (Contents*lbtoslug)/Volume; // slugs/in^3 + + } + char property_name[80]; snprintf(property_name, 80, "propulsion/tank[%d]/contents-lbs", TankNumber); PropertyManager->Tie( property_name, (FGTank*)this, &FGTank::GetContents, @@ -160,6 +199,9 @@ double FGTank::Drain(double used) PctFull = 0.0; Selected = false; } + + if (grainType != gtUNKNOWN) CalculateInertias(); + return remaining; } @@ -210,6 +252,35 @@ double FGTank::Calculate(double dt) return Temperature += (dTemp + dTemp); // For now, assume upper/lower the same } +//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +// This function calculates the moments of inertia for a solid propellant +// grain - either an end burning cylindrical grain or a bored cylindrical +// grain. + +void FGTank::CalculateInertias(void) +{ + double Mass = Contents*lbtoslug; + double RadSumSqr; + double Rad2 = Radius*Radius; + Volume = (Contents*lbtoslug)/Density; // in^3 + + switch (grainType) { + case gtCYLINDRICAL: + InnerRadius = sqrt(Rad2 - Volume/(M_PI * Length)); + RadSumSqr = (Rad2 + InnerRadius*InnerRadius)/144.0; + Ixx = 0.5*Mass*RadSumSqr; + Iyy = Mass*(3.0*RadSumSqr + Length*Length/144.0)/12.0; + break; + case gtENDBURNING: + Length = Volume/(M_PI*Rad2); + Ixx = 0.5*Mass*Rad2/144.0; + Iyy = Mass*(3.0*Rad2 + Length*Length)/(144.0*12.0); + break; + } + Izz = Iyy; + +} + //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% // The bitmasked value choices are as follows: // unset: In this case (the default) JSBSim would only print diff --git a/src/FDM/JSBSim/models/propulsion/FGTank.h b/src/FDM/JSBSim/models/propulsion/FGTank.h index f9bb0332f..950c0c12a 100644 --- a/src/FDM/JSBSim/models/propulsion/FGTank.h +++ b/src/FDM/JSBSim/models/propulsion/FGTank.h @@ -109,6 +109,9 @@ CLASS DOCUMENTATION @code + + {number} + {number} {number} @@ -119,7 +122,7 @@ CLASS DOCUMENTATION {number} {number} - {number} + {number} {number} {number} {number} @@ -131,6 +134,8 @@ CLASS DOCUMENTATION - \b type - One of FUEL or OXIDIZER. This is required. - \b radius - Equivalent radius of tank for modeling slosh, defaults to inches. +- \b grain_config type - One of CYLINDRICAL or ENDBURNING. +- \b length - length of tank for modeling solid fuel propellant grain, defaults to inches. - \b capacity - Capacity, defaults to pounds. - \b contents - Initial contents, defaults to pounds. - \b temperature - Initial temperature, defaults to degrees Fahrenheit. @@ -236,6 +241,10 @@ public: is given, otherwise 32 degrees F is returned. */ double GetTemperature(void) {return CelsiusToFahrenheit(Temperature);} + double GetIxx(void) {return Ixx;} + double GetIyy(void) {return Iyy;} + double GetIzz(void) {return Izz;} + double GetStandpipe(void) {return Standpipe;} const FGColumnVector3 GetXYZ(void); @@ -247,15 +256,25 @@ public: void SetStandpipe(double amount) { Standpipe = amount; } enum TankType {ttUNKNOWN, ttFUEL, ttOXIDIZER}; + enum GrainType {gtUNKNOWN, gtCYLINDRICAL, gtENDBURNING}; private: TankType Type; + GrainType grainType; int TankNumber; string type; + string strGType; FGColumnVector3 vXYZ; FGColumnVector3 vXYZ_drain; double Capacity; double Radius; + double InnerRadius; + double Length; + double Volume; + double Density; + double Ixx; + double Iyy; + double Izz; double PctFull; double Contents, InitialContents; double Area; @@ -264,6 +283,7 @@ private: bool Selected; FGAuxiliary* Auxiliary; FGPropertyManager* PropertyManager; + void CalculateInertias(void); void Debug(int from); }; } diff --git a/src/FDM/JSBSim/models/propulsion/FGTurbine.cpp b/src/FDM/JSBSim/models/propulsion/FGTurbine.cpp index be8dc8a19..022d39235 100644 --- a/src/FDM/JSBSim/models/propulsion/FGTurbine.cpp +++ b/src/FDM/JSBSim/models/propulsion/FGTurbine.cpp @@ -361,7 +361,10 @@ double FGTurbine::Trim() double FGTurbine::CalcFuelNeed(void) { - return FuelFlow_pph /3600 * State->Getdt() * Propulsion->GetRate(); + double dT = State->Getdt() * Propulsion->GetRate(); + FuelFlowRate = FuelFlow_pph / 3600.0; // Calculates flow in lbs/sec from lbs/hr + FuelExpended = FuelFlowRate * dT; // Calculates fuel expended in this time step + return FuelExpended; } //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% @@ -493,8 +496,6 @@ void FGTurbine::bindmodel() PropertyManager->Tie( property_name, &N1); snprintf(property_name, 80, "propulsion/engine[%u]/n2", EngineNumber); PropertyManager->Tie( property_name, &N2); - snprintf(property_name, 80, "propulsion/engine[%u]/thrust", EngineNumber); - PropertyManager->Tie( property_name, this, &FGTurbine::GetThrust); snprintf(property_name, 80, "propulsion/engine[%u]/injection_cmd", EngineNumber); PropertyManager->Tie( property_name, (FGTurbine*)this, &FGTurbine::GetInjection, &FGTurbine::SetInjection); diff --git a/src/FDM/JSBSim/models/propulsion/FGTurbine.h b/src/FDM/JSBSim/models/propulsion/FGTurbine.h index 752c622c3..52193db99 100644 --- a/src/FDM/JSBSim/models/propulsion/FGTurbine.h +++ b/src/FDM/JSBSim/models/propulsion/FGTurbine.h @@ -169,7 +169,7 @@ public: double Calculate(void); double CalcFuelNeed(void); double GetPowerAvailable(void); - double GetThrust(void) const {return Thrust;} +// double GetThrust(void) const {return Thrust;} /** A lag filter. Used to control the rate at which values are allowed to change. @param var a pointer to a variable of type double diff --git a/src/FDM/JSBSim/models/propulsion/FGTurboProp.cpp b/src/FDM/JSBSim/models/propulsion/FGTurboProp.cpp index d7be93c91..be304747f 100755 --- a/src/FDM/JSBSim/models/propulsion/FGTurboProp.cpp +++ b/src/FDM/JSBSim/models/propulsion/FGTurboProp.cpp @@ -398,7 +398,10 @@ double FGTurboProp::Start(void) double FGTurboProp::CalcFuelNeed(void) { - return FuelFlow_pph /3600 * State->Getdt() * Propulsion->GetRate(); + double dT = State->Getdt() * Propulsion->GetRate(); + FuelFlowRate = FuelFlow_pph / 3600.0; + FuelExpended = FuelFlowRate * dT; + return FuelExpended; } //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -- 2.39.5