1
0
Fork 0

Sync. with JSBSim CVS

This commit is contained in:
ehofman 2008-11-30 10:44:29 +00:00
parent 98b3701655
commit b7ebc7d78d
29 changed files with 640 additions and 610 deletions

View file

@ -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)

View file

@ -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

View file

@ -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"

View file

@ -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$

View file

@ -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) {

View file

@ -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.
<i>This function may <b>only</b> be used if Run() is called first.</i> */
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 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; }
/** 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;
// 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.
inline double GetGustNED(int idx) const {return vGustNED(idx);}
/// Retrieves the gust components in NED frame.
inline FGColumnVector3& GetGustNED(void) {return vGustNED;}
/** 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; }
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);

View file

@ -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));
}

View file

@ -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();

View file

@ -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.back()->bind(PropertyManager, PointMasses.size()-1);
}
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
double FGMassBalance::GetPointMassWeight(void)
double FGMassBalance::GetTotalPointMassWeight(void)
{
double PM_total_weight = 0.0;

View file

@ -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 <struct PointMass*> PointMasses;

View file

@ -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));

View file

@ -462,9 +462,13 @@ FGMatrix33& FGPropulsion::CalculateTankInertias(void)
tankJ = FGMatrix33();
for (unsigned int i=0; i<size; i++)
for (unsigned int i=0; i<size; i++) {
tankJ += MassBalance->GetPointmassInertia( 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;
}

View file

@ -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) {

View file

@ -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) {

View file

@ -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.
<i>This function may <b>only</b> be used if Run() is called first.</i> */
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;
void Calculate(double altitude);
void Turbulence(void);
void Debug(int from);
};

View file

@ -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;

View file

@ -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; i<SourceTanks.size(); i++) {
Tank = Propulsion->GetTank(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;
}
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

View file

@ -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;

View file

@ -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

View file

@ -63,18 +63,14 @@ CLASS DOCUMENTATION
@code
<nozzle name="{string}">
<pe unit="{PSF}"> {number} </pe>
<expr> {number} </expr>
<nzl_eff> {number} </nzl_eff>
<diam unit="{FT | M | IN}"> {number} </diam>
<area unit="{FT2 | M2 | IN2}"> {number} </area>
</nozzle>
@endcode
<h3>Configuration parameters are:</h3>
<pre>
<b>pe</b> - Nozzle exit pressure.
<b>expr</b> - Nozzle expansion ratio, Ae/At, sqft. dimensionless ratio.
<b>nzl_eff</b> - Nozzle efficiency, 0.0 - 1.0.
<b>diam</b> - Nozzle diameter.
<b>pe</b> - Nozzle design exit pressure.
<b>area</b> - Nozzle area at the exit plane.
</pre>
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);
};
}

View file

@ -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
*/

View file

@ -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
<piston_engine name="{string}">
<minmp unit="{INHG | PA | ATM}"> {number} </minmp>
<maxmp unit="{INHG | PA | ATM}"> {number} </maxmp>
<minmp unit="{INHG | PA | ATM}"> {number} </minmp> <!-- Depricated -->
<maxmp unit="{INHG | PA | ATM}"> {number} </maxmp> <!-- Depricated -->
<displacement unit="{IN3 | LTR | CC}"> {number} </displacement>
<sparkfaildrop> {number} </sparkfaildrop>
<maxhp unit="{HP | WATTS}"> {number} </maxhp>
<cycles> {number} </cycles>
<idlerpm> {number} </idlerpm>
<maxrpm> {number} </maxrpm>
<maxthrottle> {number} </maxthrottle>
<minthrottle> {number} </minthrottle>
<numboostspeeds> {number} </numboostspeeds>
<bsfc unit="{LBS/HP*HR | "KG/KW*HR"}"> {number} </bsft>
<volumetric_efficiency> {number} </volumetric_efficiency>
<boostoverride> {0 | 1} </boostoverride>
<ratedboost1 unit="{INHG | PA | ATM}"> {number} </ratedboost1>
<ratedpower1 unit="{HP | WATTS}"> {number} </ratedpower1>
@ -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).

View file

@ -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) {
if (Throttle == 1 || BurnTime > 0.0) {
if (ThrustTable != 0L) { // Thrust table given -> Solid fuel used
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; i<SourceTanks.size(); i++) {
Tank = Propulsion->GetTank(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;
}
Throttle = ThrustTable->GetValue(BurnTime);
}
if (Throttle < MinThrottle || Starved) {
PctPower = Thrust = 0.0; // desired thrust
Flameout = true;
PC = 0.0;
// 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; i<SourceTanks.size(); i++) {
Tank = Propulsion->GetTank(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;
}
}
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

View file

@ -61,10 +61,7 @@ CLASS DOCUMENTATION
/** Models a generic rocket engine.
The rocket engine is modeled given the following parameters:
<ul>
<li>Chamber pressure (in psf)</li>
<li>Specific heat ratio (usually about 1.2 for hydrocarbon fuel and LOX)</li>
<li>Propulsive efficiency (in percent, from 0 to 1.0)</li>
<li>Variance (in percent, from 0 to 1.0, nominally 0.05)</li>
<li>Specific Impulse (in sec)</li>
</ul>
Additionally, the following control inputs, operating characteristics, and
location are required, as with all other engine types:
@ -78,12 +75,8 @@ CLASS DOCUMENTATION
<li>Pitch and Yaw</li>
</ul>
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;

View file

@ -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

View file

@ -109,6 +109,9 @@ CLASS DOCUMENTATION
@code
<tank type="{FUEL | OXIDIZER}">
<grain_config type="{CYLINDRICAL | ENDBURNING}">
<length unit="{IN | FT | M}"> {number} </radius>
</grain_config>
<location unit="{FT | M | IN}">
<x> {number} </x>
<y> {number} </y>
@ -119,7 +122,7 @@ CLASS DOCUMENTATION
<y> {number} </y>
<z> {number} </z>
</drain_location>
<radius unit="{FT | M}"> {number} </radius>
<radius unit="{IN | FT | M}"> {number} </radius>
<capacity unit="{LBS | KG}"> {number} </capacity>
<contents unit="{LBS | KG}"> {number} </contents>
<temperature> {number} </temperature> <!-- must be degrees fahrenheit -->
@ -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);
};
}

View file

@ -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);

View file

@ -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

View file

@ -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;
}
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%