models/propulsion/FGRotor.h
models/propulsion/FGTank.h
models/propulsion/FGThruster.h
+ models/propulsion/FGTransmission.h
models/propulsion/FGTurbine.h
models/propulsion/FGTurboProp.h
)
models/propulsion/FGRotor.cpp
models/propulsion/FGTank.cpp
models/propulsion/FGThruster.cpp
+ models/propulsion/FGTransmission.cpp
models/propulsion/FGTurbine.cpp
models/propulsion/FGTurboProp.cpp
)
namespace JSBSim {
-static const char *IdSrc = "$Id: FGFDMExec.cpp,v 1.120 2011/11/10 12:06:13 jberndt Exp $";
+static const char *IdSrc = "$Id: FGFDMExec.cpp,v 1.133 2012/04/14 18:10:43 bcoconni Exp $";
static const char *IdHdr = ID_FDMEXEC;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
StandAlone = false;
firstPass = true;
+ IncrementThenHolding = false; // increment then hold is off by default
+ TimeStepsUntilHold = -1;
+
sim_time = 0.0;
dT = 1.0/120.0; // a default timestep size. This is needed for when JSBSim is
// run in standalone mode with no initialization file.
// Store this FDM's ID
IdFDM = (*FDMctr); // The main (parent) JSBSim instance is always the "zeroth"
-
+
// Prepare FDMctr for the next child FDM id
(*FDMctr)++; // instance. "child" instances are loaded last.
Constructing = true;
typedef int (FGFDMExec::*iPMF)(void) const;
+// typedef unsigned int (FGFDMExec::*uiPMF)(void) const;
// instance->Tie("simulation/do_trim_analysis", this, (iPMF)0, &FGFDMExec::DoTrimAnalysis, false);
instance->Tie("simulation/do_simple_trim", this, (iPMF)0, &FGFDMExec::DoTrim, false);
instance->Tie("simulation/reset", this, (iPMF)0, &FGFDMExec::ResetToInitialConditions, false);
+ instance->Tie("simulation/randomseed", this, (iPMF)0, &FGFDMExec::SRand, false);
instance->Tie("simulation/terminate", (int *)&Terminate);
instance->Tie("simulation/sim-time-sec", this, &FGFDMExec::GetSimTime);
instance->Tie("simulation/jsbsim-debug", this, &FGFDMExec::GetDebugLevel, &FGFDMExec::SetDebugLevel);
try {
Unbind();
DeAllocate();
-
+
if (IdFDM == 0) { // Meaning this is no child FDM
if(Root != 0) {
if(StandAlone)
Propulsion->in.PropFeather = FCS->GetPropFeather();
Propulsion->in.H_agl = Propagate->GetDistanceAGL();
Propulsion->in.PQR = Propagate->GetPQR();
+ Propulsion->in.vXYZcg = MassBalance->GetXYZcg();
+
break;
case eAerodynamics:
Aerodynamics->in.Alpha = Auxiliary->Getalpha();
GroundReactions->in.TotalDeltaT = dT * GroundReactions->GetRate();
GroundReactions->in.WOW = GroundReactions->GetWOW();
GroundReactions->in.Location = Propagate->GetLocation();
- for (int i=0; i<GroundReactions->GetNumGearUnits(); i++) {
- GroundReactions->in.vWhlBodyVec[i] = MassBalance->StructuralToBody(GroundReactions->GetGearUnit(i)->GetLocation());
- }
+ GroundReactions->in.vXYZcg = MassBalance->GetXYZcg();
break;
case eExternalReactions:
// There are no external inputs to this model.
Accelerations->in.Ti2b = Propagate->GetTi2b();
Accelerations->in.Tb2i = Propagate->GetTb2i();
Accelerations->in.Tec2b = Propagate->GetTec2b();
- Accelerations->in.Tl2b = Propagate->GetTl2b();
+ Accelerations->in.Tec2i = Propagate->GetTec2i();
Accelerations->in.qAttitudeECI = Propagate->GetQuaternionECI();
Accelerations->in.Moment = Aircraft->GetMoments();
Accelerations->in.GroundMoment = GroundReactions->GetMoments();
Aerodynamics->in.Wingspan = Aircraft->GetWingSpan();
Auxiliary->in.Wingspan = Aircraft->GetWingSpan();
Auxiliary->in.Wingchord = Aircraft->Getcbar();
- for (int i=0; i<GroundReactions->GetNumGearUnits(); i++) {
- GroundReactions->in.vWhlBodyVec[i] = MassBalance->StructuralToBody(GroundReactions->GetGearUnit(i)->GetLocation());
- }
+ GroundReactions->in.vXYZcg = MassBalance->GetXYZcg();
+ Propulsion->in.vXYZcg = MassBalance->GetXYZcg();
LoadPlanetConstants();
}
{
if (mode == 1) {
for (unsigned int i=0; i<Outputs.size(); i++) {
- Outputs[i]->SetStartNewFile(true);
+ Outputs[i]->SetStartNewFile(true);
}
}
-
+
ResetToInitialConditions();
}
while (element) {
if (debug_lvl > 0) cout << endl << " Output data set: " << idx << " ";
FGOutput* Output = new FGOutput(this);
- Output->InitModel();
- Schedule(Output);
result = Output->Load(element);
if (!result) {
cerr << endl << "Aircraft output element has problems in file " << aircraftCfgFileName << endl;
+ delete Output;
return result;
} else {
+ Output->InitModel();
+ Schedule(Output);
Outputs.push_back(Output);
string outputProp = CreateIndexedPropertyName("simulation/output",idx);
instance->Tie(outputProp+"/log_rate_hz", Output, (iOPMF)0, &FGOutput::SetRate, false);
- instance->Tie("simulation/force-output", this, (iOPV)0, &FGFDMExec::ForceOutput, false);
idx++;
}
element = document->FindNextElement("output");
}
+ if (idx)
+ instance->Tie("simulation/force-output", this, (iOPV)0, &FGFDMExec::ForceOutput, false);
// Lastly, process the child element. This element is OPTIONAL - and NOT YET SUPPORTED.
element = document->FindElement("child");
<< "-------------------------------------------------------------------------------"
<< reset << endl;
}
-
+
if (IsChild) debug_lvl = saved_debug_lvl;
} else {
cerr << endl << highint << fgred << " No location was found for this child object!" << reset << endl;
exit(-1);
}
-
+
Element* orientation = el->FindElement("orient");
if (orientation) {
child->Orient = orientation->FindElementTripletConvertTo("RAD");
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+void FGFDMExec::CheckIncrementalHold(void)
+{
+ // Only check if increment then hold is on
+ if( IncrementThenHolding ) {
+
+ if (TimeStepsUntilHold == 0) {
+
+ // Should hold simulation if TimeStepsUntilHold has reached zero
+ holding = true;
+
+ // Still need to decrement TimeStepsUntilHold as value of -1
+ // indicates that incremental then hold is turned off
+ IncrementThenHolding = false;
+ TimeStepsUntilHold--;
+
+ } else if ( TimeStepsUntilHold > 0 ) {
+ // Keep decrementing until 0 is reached
+ TimeStepsUntilHold--;
+ }
+ }
+}
+
+//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
void FGFDMExec::ForceOutput(int idx)
{
if (idx >= (int)0 && idx < (int)Outputs.size()) Outputs[idx]->Print();
FGOutput* Output = new FGOutput(this);
Output->SetDirectivesFile(RootDir + fname);
- Output->InitModel();
- Schedule(Output);
result = Output->Load(0);
if (result) {
+ Output->InitModel();
+ Schedule(Output);
Output->Run(holding);
Outputs.push_back(Output);
typedef double (FGOutput::*iOPMF)(void) const;
sim_time = saved_time;
}
+//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+void FGFDMExec::SRand(int sr)
+{
+ srand(sr);
+}
+
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
// The bitmasked value choices are as follows:
// unset: In this case (the default) JSBSim would only print
if (debug_lvl & 1 && IdFDM == 0) { // Standard console startup message output
if (from == 0) { // Constructor
- cout << "\n\n "
+ cout << "\n\n "
<< "JSBSim Flight Dynamics Model v" << JSBSim_version << endl;
cout << " [JSBSim-ML v" << needed_cfg_version << "]\n\n";
cout << "JSBSim startup beginning ...\n\n";
DEFINITIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#define ID_FDMEXEC "$Id: FGFDMExec.h,v 1.74 2011/11/09 21:58:26 bcoconni Exp $"
+#define ID_FDMEXEC "$Id: FGFDMExec.h,v 1.76 2012/01/21 16:46:08 jberndt Exp $"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FORWARD DECLARATIONS
property actually maps toa function call of DoTrim().
@author Jon S. Berndt
- @version $Revision: 1.74 $
+ @version $Revision: 1.76 $
*/
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
void EnableOutput(void);
/// Pauses execution by preventing time from incrementing.
void Hold(void) {holding = true;}
+ /// Turn on hold after increment
+ void EnableIncrementThenHold(int Timesteps) {TimeStepsUntilHold = Timesteps; IncrementThenHolding = true;}
+ /// Checks if required to hold afer increment
+ void CheckIncrementalHold(void);
/// Resumes execution from a "Hold".
void Resume(void) {holding = false;}
/// Returns true if the simulation is Holding (i.e. simulation time is not moving).
double saved_dT;
double sim_time;
bool holding;
+ bool IncrementThenHolding;
+ int TimeStepsUntilHold;
bool Constructing;
bool modelLoaded;
bool IsChild;
bool ReadChild(Element*);
bool ReadPrologue(Element*);
void ResetToInitialConditions(int mode);
+ void SRand(int sr);
void LoadInputs(unsigned int idx);
void LoadPlanetConstants(void);
void LoadModelConstants(void);
namespace JSBSim {
-static const char *IdSrc = "$Id: FGJSBBase.cpp,v 1.32 2011/10/22 14:38:30 bcoconni Exp $";
+static const char *IdSrc = "$Id: FGJSBBase.cpp,v 1.35 2012/03/25 11:05:36 bcoconni Exp $";
static const char *IdHdr = ID_JSBBASE;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
msg.messageId = messageId++;
msg.subsystem = "FDM";
msg.type = Message::eInteger;
- msg.bVal = (iVal != 0);
+ msg.iVal = iVal;
Messages.push(msg);
}
msg.messageId = messageId++;
msg.subsystem = "FDM";
msg.type = Message::eDouble;
- msg.bVal = (dVal != 0.0);
+ msg.dVal = dVal;
Messages.push(msg);
}
if (fgGetBool("/environment/params/control-fdm-atmosphere")) {
Atmosphere->SetTemperature(temperature->getDoubleValue(), get_Altitude(), FGAtmosphere::eCelsius);
- Atmosphere->SetPressureSL(pressureSL->getDoubleValue(), FGAtmosphere::eInchesHg);
+ Atmosphere->SetPressureSL(FGAtmosphere::eInchesHg, pressureSL->getDoubleValue());
// initialize to no turbulence, these values get set in the update loop
Winds->SetTurbType(FGWinds::ttNone);
Winds->SetTurbGain(0.0);
}
Atmosphere->SetTemperature(temperature->getDoubleValue(), get_Altitude(), FGAtmosphere::eCelsius);
- Atmosphere->SetPressureSL(pressureSL->getDoubleValue(), FGAtmosphere::eInchesHg);
+ Atmosphere->SetPressureSL(FGAtmosphere::eInchesHg, pressureSL->getDoubleValue());
Winds->SetTurbType((FGWinds::tType)TURBULENCE_TYPE_NAMES[turbulence_model->getStringValue()]);
switch( Winds->GetTurbType() ) {
node->setDoubleValue("rollspeed-ms", gear->GetWheelRollVel()*0.3043);
node->setBoolValue("has-brake", gear->GetBrakeGroup() > 0);
node->setDoubleValue("position-norm", gear->GetGearUnitPos());
- node->setDoubleValue("tire-pressure-norm", gear->GetTirePressure());
+// node->setDoubleValue("tire-pressure-norm", gear->GetTirePressure());
node->setDoubleValue("compression-norm", gear->GetCompLen());
node->setDoubleValue("compression-ft", gear->GetCompLen());
if ( gear->GetSteerable() )
node->getChild("wow", 0, true)->setBoolValue( gear->GetWOW());
node->getChild("rollspeed-ms", 0, true)->setDoubleValue(gear->GetWheelRollVel()*0.3043);
node->getChild("position-norm", 0, true)->setDoubleValue(gear->GetGearUnitPos());
- gear->SetTirePressure(node->getDoubleValue("tire-pressure-norm"));
+// gear->SetTirePressure(node->getDoubleValue("tire-pressure-norm"));
node->setDoubleValue("compression-norm", gear->GetCompLen());
node->setDoubleValue("compression-ft", gear->GetCompLen());
if ( gear->GetSteerable() )
namespace JSBSim {
-static const char *IdSrc = "$Id: FGInitialCondition.cpp,v 1.78 2011/11/09 21:57:51 bcoconni Exp $";
+static const char *IdSrc = "$Id: FGInitialCondition.cpp,v 1.81 2012/04/08 15:22:56 jberndt Exp $";
static const char *IdHdr = ID_INITIALCONDITION;
//******************************************************************************
filename = "initfile.xml";
else
filename.append("/initfile.xml");
-
+
ofstream outfile(filename.c_str());
FGPropagate* Propagate = fdmex->GetPropagate();
-
+
if (outfile.is_open()) {
outfile << "<?xml version=\"1.0\"?>" << endl;
outfile << "<initialize name=\"reset00\">" << endl;
// Refer to Stevens and Lewis, 1.5-14a, pg. 49.
// This is the rotation rate of the "Local" frame, expressed in the local frame.
+ const FGMatrix33& Tl2b = orientation.GetT();
double radInv = 1.0 / position.GetRadius();
FGColumnVector3 vOmegaLocal = FGColumnVector3(
radInv*vUVW_NED(eEast),
-radInv*vUVW_NED(eNorth),
-radInv*vUVW_NED(eEast)*position.GetTanLatitude() );
- vPQR_body = vOmegaLocal;
+ vPQR_body = Tl2b * vOmegaLocal;
return result;
}
// - Local
// - Body
// The vehicle will be defaulted to (0,0,0) in the Body frame if nothing is provided.
-
+
Element* velocity_el = document->FindElement("velocity");
FGColumnVector3 vInitVelocity = FGColumnVector3(0.0, 0.0, 0.0);
FGMatrix33 mTec2l = position.GetTec2l();
// - ECI (Earth Centered Inertial)
// - ECEF (Earth Centered, Earth Fixed)
// - Body
-
+
Element* attrate_el = document->FindElement("attitude_rate");
const FGMatrix33& Tl2b = orientation.GetT();
} else if (frame == "ecef") {
vPQR_body = Tl2b * position.GetTec2l() * vAttRate;
} else if (frame == "local") {
- vPQR_body = vAttRate + vOmegaLocal;
+ vPQR_body = Tl2b * (vAttRate + vOmegaLocal);
+ } else if (frame == "body") {
+ vPQR_body = vAttRate;
} else if (!frame.empty()) { // misspelling of frame
-
+
cerr << endl << fgred << " Attitude rate frame type: \"" << frame
<< "\" is not supported!" << reset << endl << endl;
result = false;
} else if (frame.empty()) {
- vPQR_body = vOmegaLocal;
+ vPQR_body = Tl2b * vOmegaLocal;
}
} else { // Body frame attitude rate assumed 0 relative to local.
- vPQR_body = vOmegaLocal;
+ vPQR_body = Tl2b * vOmegaLocal;
}
return result;
#include "input_output/FGXMLFileRead.h"
#include "math/FGLocation.h"
+#include "math/FGQuaternion.h"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
DEFINITIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#define ID_INITIALCONDITION "$Id: FGInitialCondition.h,v 1.32 2011/11/06 18:14:51 bcoconni Exp $"
+#define ID_INITIALCONDITION "$Id: FGInitialCondition.h,v 1.33 2011/12/10 15:49:21 bcoconni Exp $"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FORWARD DECLARATIONS
@property ic/r-rad_sec (read/write) Yaw rate initial condition in radians/second
@author Tony Peden
- @version "$Id: FGInitialCondition.h,v 1.32 2011/11/06 18:14:51 bcoconni Exp $"
+ @version "$Id: FGInitialCondition.h,v 1.33 2011/12/10 15:49:21 bcoconni Exp $"
*/
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
namespace JSBSim {
-static const char *IdSrc = "$Id: FGColumnVector3.cpp,v 1.14 2010/12/07 12:57:14 jberndt Exp $";
+static const char *IdSrc = "$Id: FGColumnVector3.cpp,v 1.15 2012/02/07 00:27:51 jentron Exp $";
static const char *IdHdr = ID_COLUMNVECTOR3;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
string FGColumnVector3::Dump(const string& delimiter) const
{
ostringstream buffer;
- buffer << std::setw(18) << std::setprecision(16) << data[0] << delimiter;
- buffer << std::setw(18) << std::setprecision(16) << data[1] << delimiter;
- buffer << std::setw(18) << std::setprecision(16) << data[2];
+ buffer << std::setprecision(16) << data[0] << delimiter;
+ buffer << std::setprecision(16) << data[1] << delimiter;
+ buffer << std::setprecision(16) << data[2];
return buffer.str();
}
DEFINITIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#define ID_COLUMNVECTOR3 "$Id: FGColumnVector3.h,v 1.16 2010/12/07 12:57:14 jberndt Exp $"
+#define ID_COLUMNVECTOR3 "$Id: FGColumnVector3.h,v 1.17 2011/12/10 15:49:21 bcoconni Exp $"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FORWARD DECLARATIONS
/** This class implements a 3 element column vector.
@author Jon S. Berndt, Tony Peden, et. al.
- @version $Id: FGColumnVector3.h,v 1.16 2010/12/07 12:57:14 jberndt Exp $
+ @version $Id: FGColumnVector3.h,v 1.17 2011/12/10 15:49:21 bcoconni Exp $
*/
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
/** Dot product of two vectors
Compute and return the euclidean dot (or scalar) product of two vectors
v1 and v2 */
- friend inline double DotProduct(const FGColumnVector3& v1, const FGColumnVector3& v2) {
+ friend double DotProduct(const FGColumnVector3& v1, const FGColumnVector3& v2) {
return v1.data[0]*v2.data[0] + v1.data[1]*v2.data[1] + v1.data[2]*v2.data[2];
}
namespace JSBSim {
-static const char *IdSrc = "$Id: FGFunction.cpp,v 1.42 2011/09/07 02:36:04 jberndt Exp $";
+static const char *IdSrc = "$Id: FGFunction.cpp,v 1.43 2012/02/05 11:15:54 bcoconni Exp $";
static const char *IdHdr = ID_FUNCTION;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
case eIfThen:
{
i = Parameters.size();
- if (i != 3) {
+ if (i == 3) {
if (GetBinary(temp) == 1) {
temp = Parameters[1]->GetValue();
} else {
namespace JSBSim {
-static const char *IdSrc = "$Id: FGLocation.cpp,v 1.26 2011/11/06 18:14:51 bcoconni Exp $";
+static const char *IdSrc = "$Id: FGLocation.cpp,v 1.29 2012/04/14 12:14:37 bcoconni Exp $";
static const char *IdHdr = ID_LOCATION;
using std::cerr;
using std::endl;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
FGLocation::FGLocation(void)
+ : mECLoc(1.0, 0.0, 0.0), mCacheValid(false)
{
- mCacheValid = false;
-
a = b = a2 = b2 = 0.0;
e = e2 = f = 1.0;
eps2 = -1.0;
mTec2i.InitMatrix();
mTi2l.InitMatrix();
mTl2i.InitMatrix();
- mECLoc.InitMatrix();
}
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FGLocation::FGLocation(double lon, double lat, double radius)
+ : mCacheValid(false)
{
- mCacheValid = false;
-
a = b = a2 = b2 = 0.0;
e = e2 = f = 1.0;
eps2 = -1.0;
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-FGLocation::FGLocation(const FGColumnVector3& lv) : mECLoc(lv), mCacheValid(false)
+FGLocation::FGLocation(const FGColumnVector3& lv)
+ : mECLoc(lv), mCacheValid(false)
{
a = b = a2 = b2 = 0.0;
e = e2 = f = 1.0;
FGLocation::FGLocation(const FGLocation& l)
: mECLoc(l.mECLoc), mCacheValid(l.mCacheValid)
{
-
a = l.a;
b = l.b;
a2 = l.a2;
void FGLocation::SetPositionGeodetic(double lon, double lat, double height)
{
mCacheValid = false;
-
+
mGeodLat = lat;
mLon = lon;
GeodeticAltitude = height;
f = 1.0 - b/a;
}
-//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-// Compute the ECEF to ECI transformation matrix using Stevens and Lewis "Aircraft
-// Control and Simulation", second edition, eqn. 1.4-12, pg. 39. In Stevens and Lewis
-// notation, this is C_i/e, a transformation from ECEF to ECI.
-
-const FGMatrix33& FGLocation::GetTec2i(void)
-{
- ComputeDerived();
- return mTec2i;
-}
-
-//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-// This is given in Stevens and Lewis "Aircraft
-// Control and Simulation", second edition, eqn. 1.4-12, pg. 39
-// The notation in Stevens and Lewis is: C_e/i. This represents a transformation
-// from ECI to ECEF - and the orientation of the ECEF frame relative to the ECI
-// frame.
-
-const FGMatrix33& FGLocation::GetTi2ec(void)
-{
- ComputeDerived();
- return mTi2ec;
-}
-
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
void FGLocation::ComputeDerivedUnconditional(void) const
// Compute the transform matrices from and to the earth centered frame.
// See Stevens and Lewis, "Aircraft Control and Simulation", Second Edition,
- // Eqn. 1.4-13, page 40. In Stevens and Lewis notation, this is C_n/e - the
+ // Eqn. 1.4-13, page 40. In Stevens and Lewis notation, this is C_n/e - the
// orientation of the navigation (local) frame relative to the ECEF frame,
// and a transformation from ECEF to nav (local) frame.
-sinLon , cosLon , 0.0 ,
-cosLon*cosLat, -sinLon*cosLat, -sinLat );
- // In Stevens and Lewis notation, this is C_e/n - the
+ // In Stevens and Lewis notation, this is C_e/n - the
// orientation of the ECEF frame relative to the nav (local) frame,
// and a transformation from nav (local) to ECEF frame.
if (a != 0.0 && b != 0.0) {
double c, p, q, s, t, u, v, w, z, p2, u2, r0;
- double Ne, P, Q0, Q, signz0, sqrt_q, z_term;
+ double Ne, P, Q0, Q, signz0, sqrt_q, z_term;
p = fabs(mECLoc(eZ))/eps2;
s = r02/(e2*eps2);
p2 = p*p;
DEFINITIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#define ID_LOCATION "$Id: FGLocation.h,v 1.29 2011/11/06 18:14:51 bcoconni Exp $"
+#define ID_LOCATION "$Id: FGLocation.h,v 1.31 2012/02/05 14:56:17 bcoconni Exp $"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FORWARD DECLARATIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
/** FGLocation holds an arbitrary location in the Earth centered Earth fixed
- reference frame (ECEF). This coordinate frame has its center in the middle
+ reference frame (ECEF). The coordinate frame ECEF has its center in the middle
of the earth. The X-axis points from the center of the Earth towards a
location with zero latitude and longitude on the Earth surface. The Y-axis
points from the center of the Earth towards a location with zero latitude
conversion functions for conversion of position vectors given in the one
frame to positions in the other frame.
+ To keep the transformation matrices between the ECI and ECEF frames up to
+ date, the Earth angular position must be updated by calling
+ SetEarthPositionAngle() or IncrementEarthPositionAngle(). This must be done
+ prior to any conversion from and to the ECI frame.
+
The Earth centered reference frame is NOT an inertial frame since it rotates
with the Earth.
- The coordinates in the Earth centered frame are the master values. All other
+ The cartesian coordinates (X,Y,Z) in the Earth centered frame are the master values. All other
values are computed from these master values and are cached as long as the
location is changed by access through a non-const member function. Values
are cached to improve performance. It is best practice to work with a
Given,
- -that we model a vehicle near the Earth
- -that the Earth surface radius is about 2*10^7, ft
- -that we use double values for the representation of the location
-
+ - that we model a vehicle near the Earth
+ - that the Earth surface radius is about 2*10^7, ft
+ - that we use double values for the representation of the location
+
we have an accuracy of about
-
+
1e-16*2e7ft/1 = 2e-9 ft
-
+
left. This should be sufficient for our needs. Note that this is the same
relative accuracy we would have when we compute directly with
lon/lat/radius. For the radius value this is clear. For the lon/lat pair
@see W. C. Durham "Aircraft Dynamics & Control", section 2.2
@author Mathias Froehlich
- @version $Id: FGLocation.h,v 1.29 2011/11/06 18:14:51 bcoconni Exp $
+ @version $Id: FGLocation.h,v 1.31 2012/02/05 14:56:17 bcoconni Exp $
*/
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
@param radius distance from center of earth to vehicle in feet*/
FGLocation(double lon, double lat, double radius);
- /** Column constructor. */
+ /** Constructor to initialize the location with the cartesian coordinates
+ (X,Y,Z) contained in the input FGColumnVector3. Distances are in feet,
+ the position is expressed in the ECEF frame.
+ @param lv vector that contain the cartesian coordinates*/
FGLocation(const FGColumnVector3& lv);
/** Copy constructor. */
return -mTec2l(3,3)/cLat;
}
+ /** Return the Earth Position Angle.
+ This is the relative orientation of the ECEF frame with respect to the
+ Inertial frame.
+ @return the Earth fixed frame (ECEF) rotation offset about the axis with
+ respect to the Inertial (ECI) frame in radians. */
double GetEPA() const {return epa;}
/** Get the distance from the center of the earth.
//double GetRadius() const { return mECLoc.Magnitude(); } // may not work with FlightGear
double GetRadius() const { ComputeDerived(); return mRadius; }
- /// @name Functions that need the ground callback to be set
+ /** @name Functions that rely on the ground callback
+ The following functions allow to set and get the vehicle position above
+ the sea or the ground. The sea and the ground levels are obtained by
+ interrogating an FGGroundCallback instance. A ground callback must
+ therefore be set with SetGroundCallback() before calling any of these
+ functions. */
///@{
/** Set the altitude above sea level.
- @param altitudeASL altitude above Sea Level in feet. */
+ @param altitudeASL altitude above Sea Level in feet.
+ @see SetGroundCallback */
void SetAltitudeASL(double altitudeASL)
{ SetRadius(GroundCallback->GetSeaLevelRadius(*this) + altitudeASL); }
/** Set the altitude above ground level.
- @param altitudeAGL altitude above Ground Level in feet. */
+ @param altitudeAGL altitude above Ground Level in feet.
+ @see SetGroundCallback */
void SetAltitudeAGL(double altitudeAGL, double time)
{ SetRadius(GroundCallback->GetTerrainGeoCentRadius(time, *this) + altitudeAGL); }
/** Get the local sea level radius
- @return the sea level radius at the location in feet. */
+ @return the sea level radius at the location in feet.
+ @see SetGroundCallback */
double GetSeaLevelRadius(void) const
{ ComputeDerived(); return GroundCallback->GetSeaLevelRadius(*this); }
/** Get the local terrain radius
- @return the terrain level radius at the location in feet. */
+ @return the terrain level radius at the location in feet.
+ @see SetGroundCallback */
double GetTerrainRadius(double time) const
{ ComputeDerived(); return GroundCallback->GetTerrainGeoCentRadius(time, *this); }
/** Get the altitude above sea level.
- @return the altitude ASL in feet. */
+ @return the altitude ASL in feet.
+ @see SetGroundCallback */
double GetAltitudeASL() const
{ ComputeDerived(); return GroundCallback->GetAltitude(*this); }
/** Get the altitude above ground level.
- @return the altitude AGL in feet. */
+ @return the altitude AGL in feet.
+ @see SetGroundCallback */
double GetAltitudeAGL(double time) const {
FGLocation c;
FGColumnVector3 n,v,w;
@param normal Terrain normal vector in contact point (ECEF frame)
@param v Terrain linear velocity in contact point (ECEF frame)
@param w Terrain angular velocity in contact point (ECEF frame)
- @return Location altitude above contact point (AGL) in feet. */
+ @return Location altitude above contact point (AGL) in feet.
+ @see SetGroundCallback */
double GetContactPoint(double time,
FGLocation& contact, FGColumnVector3& normal,
FGColumnVector3& v, FGColumnVector3& w) const
{ ComputeDerived(); return GroundCallback->GetAGLevel(time, *this, contact, normal, v, w); }
+ ///@}
- /** Sets the ground callback pointer. For optimal memory management, a shared
- pointer is used internally that maintains a reference counter. The calling
- application must therefore use FGGroundCallback_ptr 'smart pointers' to
- manage their copy of the ground callback.
+ /** Sets the ground callback pointer. The FGGroundCallback instance will be
+ interrogated by FGLocation each time some terrain informations are needed.
+ This will mainly occur when altitudes above the sea level or above the
+ ground level are needed. A 'smart pointer' is used internally to prevent
+ the FGGroundCallback instance against accidental deletion. This can only
+ work if the calling application also make use of FGGroundCallback_ptr
+ 'smart pointers' to manage their copy of the ground callback.
@param gc A pointer to a ground callback object
@see FGGroundCallback
*/
static void SetGroundCallback(FGGroundCallback* gc) { GroundCallback = gc; }
- /** Get a pointer to the ground callback currently used. It is recommanded
- to store the returned pointer in a 'smart pointer' FGGroundCallback_ptr.
+ /** Get a pointer to the ground callback currently used. Since the
+ FGGroundcallback instance might have been created outside JSBSim, it is
+ recommanded to store the returned pointer in a 'smart pointer'
+ FGGroundCallback_ptr. This pointer maintains a reference counter and
+ protects the returned pointer against an accidental deletion of the object
+ it is pointing to.
@return A pointer to the current ground callback object.
@see FGGroundCallback
*/
static FGGroundCallback* GetGroundCallback(void) { return GroundCallback; }
- ///@}
/** Transform matrix from local horizontal to earth centered frame.
- Returns a const reference to the rotation matrix of the transform from
+ @return a const reference to the rotation matrix of the transform from
the local horizontal frame to the earth centered frame. */
const FGMatrix33& GetTl2ec(void) const { ComputeDerived(); return mTl2ec; }
/** Transform matrix from the earth centered to local horizontal frame.
- Returns a const reference to the rotation matrix of the transform from
+ @return a const reference to the rotation matrix of the transform from
the earth centered frame to the local horizontal frame. */
const FGMatrix33& GetTec2l(void) const { ComputeDerived(); return mTec2l; }
/** Transform matrix from inertial to earth centered frame.
- Returns a const reference to the rotation matrix of the transform from
- the inertial frame to the earth centered frame (ECI to ECEF). */
- const FGMatrix33& GetTi2ec(void);
+ @return a const reference to the rotation matrix of the transform from
+ the inertial frame to the earth centered frame (ECI to ECEF).
+ @see SetEarthPositionAngle
+ @see IncrementEarthPositionAngle */
+ const FGMatrix33& GetTi2ec(void) const { ComputeDerived(); return mTi2ec; }
/** Transform matrix from the earth centered to inertial frame.
- Returns a const reference to the rotation matrix of the transform from
- the earth centered frame to the inertial frame (ECEF to ECI). */
- const FGMatrix33& GetTec2i(void);
-
+ @return a const reference to the rotation matrix of the transform from
+ the earth centered frame to the inertial frame (ECEF to ECI).
+ @see SetEarthPositionAngle
+ @see IncrementEarthPositionAngle */
+ const FGMatrix33& GetTec2i(void) const { ComputeDerived(); return mTec2i; }
+
+ /** Transform matrix from the inertial to local horizontal frame.
+ @return a const reference to the rotation matrix of the transform from
+ the inertial frame to the local horizontal frame.
+ @see SetEarthPositionAngle
+ @see IncrementEarthPositionAngle */
const FGMatrix33& GetTi2l(void) const {ComputeDerived(); return mTi2l;}
+ /** Transform matrix from local horizontal to inertial frame.
+ @return a const reference to the rotation matrix of the transform from
+ the local horizontal frame to the inertial frame.
+ @see SetEarthPositionAngle
+ @see IncrementEarthPositionAngle */
const FGMatrix33& GetTl2i(void) const {ComputeDerived(); return mTl2i;}
/** Conversion from Local frame coordinates to a location in the
earth centered and fixed frame.
+ This function calculates the FGLocation of an object which position
+ relative to the vehicle is given as in input.
@param lvec Vector in the local horizontal coordinate frame
@return The location in the earth centered and fixed frame */
FGLocation LocalToLocation(const FGColumnVector3& lvec) const {
/** Conversion from a location in the earth centered and fixed frame
to local horizontal frame coordinates.
+ This function calculates the relative position between the vehicle and
+ the input vector and returns the result expressed in the local frame.
@param ecvec Vector in the earth centered and fixed frame
@return The vector in the local horizontal coordinate frame */
FGColumnVector3 LocationToLocal(const FGColumnVector3& ecvec) const {
The location can be set by an Earth-centered, Earth-fixed (ECEF) frame
position vector. The cache is marked as invalid, so any future requests
for selected important data will cause the parameters to be calculated.
- @param v the ECEF column vector in feet.
+ @param v the ECEF column vector in feet.
@return a reference to the FGLocation object. */
const FGLocation& operator=(const FGColumnVector3& v)
{
mECLoc(eY) = v(eY);
mECLoc(eZ) = v(eZ);
mCacheValid = false;
- ComputeDerived();
+ //ComputeDerived();
return *this;
}
/** Sets this location via the supplied location object.
- @param v A location object reference.
+ @param v A location object reference.
@return a reference to the FGLocation object. */
const FGLocation& operator=(const FGLocation& l);
bool operator!=(const FGLocation& l) const { return ! operator==(l); }
/** This operator adds the ECEF position vectors.
- The supplied vector (right side) is added to the ECEF position vector
- on the left side of the equality, and a pointer to this object is
- returned. */
+ The cartesian coordinates of the supplied vector (right side) are added to
+ the ECEF position vector on the left side of the equality, and a reference
+ to this object is returned. */
const FGLocation& operator+=(const FGLocation &l) {
mCacheValid = false;
mECLoc += l.mECLoc;
return *this;
}
+ /** This operator substracts the ECEF position vectors.
+ The cartesian coordinates of the supplied vector (right side) are
+ substracted from the ECEF position vector on the left side of the
+ equality, and a reference to this object is returned. */
const FGLocation& operator-=(const FGLocation &l) {
mCacheValid = false;
mECLoc -= l.mECLoc;
return *this;
}
+ /** This operator scales the ECEF position vector.
+ The cartesian coordinates of the ECEF position vector on the left side of
+ the equality are scaled by the supplied value (right side), and a
+ reference to this object is returned. */
const FGLocation& operator*=(double scalar) {
mCacheValid = false;
mECLoc *= scalar;
return *this;
}
+ /** This operator scales the ECEF position vector.
+ The cartesian coordinates of the ECEF position vector on the left side of
+ the equality are scaled by the inverse of the supplied value (right side),
+ and a reference to this object is returned. */
const FGLocation& operator/=(double scalar) {
return operator*=(1.0/scalar);
}
+ /** This operator adds two ECEF position vectors.
+ A new object is returned that defines a position which is the sum of the
+ cartesian coordinates of the two positions provided. */
FGLocation operator+(const FGLocation& l) const {
return FGLocation(mECLoc + l.mECLoc);
}
+ /** This operator substracts two ECEF position vectors.
+ A new object is returned that defines a position which is the difference
+ of the cartesian coordinates of the two positions provided. */
FGLocation operator-(const FGLocation& l) const {
return FGLocation(mECLoc - l.mECLoc);
}
+ /** This operator scales an ECEF position vector.
+ A new object is returned that defines a position made of the cartesian
+ coordinates of the provided ECEF position scaled by the supplied scalar
+ value. */
FGLocation operator*(double scalar) const {
return FGLocation(scalar*mECLoc);
}
mutable double mRadius;
mutable double mGeodLat;
mutable double GeodeticAltitude;
-
+
double initial_longitude;
/** The cached rotation matrices from and to the associated frames. */
#include "FGMatrix33.h"
#include "FGColumnVector3.h"
+#include "FGQuaternion.h"
#include <sstream>
#include <iomanip>
namespace JSBSim {
-static const char *IdSrc = "$Id: FGMatrix33.cpp,v 1.11 2010/12/07 12:57:14 jberndt Exp $";
+static const char *IdSrc = "$Id: FGMatrix33.cpp,v 1.12 2011/12/10 15:49:21 bcoconni Exp $";
static const char *IdHdr = ID_MATRIX33;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
// Find largest of the above
idx = 0;
- for (int i=1; i<4; i++) if (tempQ[i] > tempQ[idx]) idx = i;
+ for (int i=1; i<4; i++) if (tempQ[i] > tempQ[idx]) idx = i;
switch(idx) {
case 0:
DEFINITIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#define ID_MATRIX33 "$Id: FGMatrix33.h,v 1.14 2010/12/07 12:57:14 jberndt Exp $"
+#define ID_MATRIX33 "$Id: FGMatrix33.h,v 1.15 2011/12/10 15:49:21 bcoconni Exp $"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FORWARD DECLARATIONS
/** Returns the quaternion associated with this direction cosine (rotation) matrix.
*/
FGQuaternion GetQuaternion(void);
-
+
/** Determinant of the matrix.
@return the determinant of the matrix.
*/
std::istream& operator>>(std::istream& is, FGMatrix33& M);
} // namespace JSBSim
-
-#include "FGQuaternion.h"
-
#endif
INCLUDES
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#include "FGModelFunctions.h"
+#include <iostream>
+#include <sstream>
#include <string>
+#include "FGModelFunctions.h"
using namespace std;
namespace JSBSim {
-static const char *IdSrc = "$Id: FGModelFunctions.cpp,v 1.4 2010/09/07 00:40:03 jberndt Exp $";
+static const char *IdSrc = "$Id: FGModelFunctions.cpp,v 1.5 2012/04/13 13:25:52 jberndt Exp $";
static const char *IdHdr = ID_MODELFUNCTIONS;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
(*it)->GetValue();
}
+//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+string FGModelFunctions::GetFunctionStrings(const string& delimeter) const
+{
+ string FunctionStrings = "";
+ bool firstime = true;
+ unsigned int sd;
+
+ for (sd = 0; sd < PreFunctions.size(); sd++) {
+ if (firstime) {
+ firstime = false;
+ } else {
+ FunctionStrings += delimeter;
+ }
+ FunctionStrings += PreFunctions[sd]->GetName();
+ }
+
+ for (sd = 0; sd < PostFunctions.size(); sd++) {
+ if (firstime) {
+ firstime = false;
+ } else {
+ FunctionStrings += delimeter;
+ }
+ FunctionStrings += PostFunctions[sd]->GetName();
+ }
+
+ return FunctionStrings;
+}
+
+//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+string FGModelFunctions::GetFunctionValues(const string& delimeter) const
+{
+ ostringstream buf;
+
+ for (unsigned int sd = 0; sd < PreFunctions.size(); sd++) {
+ if (buf.tellp() > 0) buf << delimeter;
+ buf << PreFunctions[sd]->GetValue();
+ }
+
+ for (unsigned int sd = 0; sd < PostFunctions.size(); sd++) {
+ if (buf.tellp() > 0) buf << delimeter;
+ buf << PostFunctions[sd]->GetValue();
+ }
+
+ return buf.str();
+}
+
+//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
}
DEFINITIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#define ID_MODELFUNCTIONS "$Id: FGModelFunctions.h,v 1.4 2011/06/21 04:41:54 jberndt Exp $"
+#define ID_MODELFUNCTIONS "$Id: FGModelFunctions.h,v 1.5 2012/04/13 13:25:52 jberndt Exp $"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FORWARD DECLARATIONS
void PreLoad(Element* el, FGPropertyManager* PropertyManager, std::string prefix="");
void PostLoad(Element* el, FGPropertyManager* PropertyManager, std::string prefix="");
+ /** Gets the strings for the current set of functions.
+ @param delimeter either a tab or comma string depending on output type
+ @return a string containing the descriptive names for all functions */
+ std::string GetFunctionStrings(const std::string& delimeter) const;
+
+ /** Gets the function values.
+ @param delimeter either a tab or comma string depending on output type
+ @return a string containing the numeric values for the current set of
+ functions */
+ std::string GetFunctionValues(const std::string& delimeter) const;
+
protected:
std::vector <FGFunction*> PreFunctions;
std::vector <FGFunction*> PostFunctions;
DEFINITIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#define ID_QUATERNION "$Id: FGQuaternion.h,v 1.23 2011/10/31 14:54:40 bcoconni Exp $"
+#define ID_QUATERNION "$Id: FGQuaternion.h,v 1.24 2011/12/10 15:49:21 bcoconni Exp $"
namespace JSBSim {
Useful for initialization of increments */
static FGQuaternion zero(void) { return FGQuaternion( 0.0, 0.0, 0.0, 0.0 ); }
+ friend FGQuaternion QExp(const FGColumnVector3& omega);
+
private:
/** Copying by assigning the vector valued components. */
FGQuaternion(double q1, double q2, double q3, double q4) : mCacheValid(false)
return FGQuaternion(scalar*q.data[0], scalar*q.data[1], scalar*q.data[2], scalar*q.data[3]);
}
+/** Quaternion exponential
+ @param omega rotation velocity
+ Calculate the unit quaternion which is the result of the exponentiation of
+ the vector 'omega'.
+*/
+inline FGQuaternion QExp(const FGColumnVector3& omega) {
+ FGQuaternion qexp;
+ double angle = omega.Magnitude();
+ double sina_a = angle > 0.0 ? sin(angle)/angle : 1.0;
+
+ qexp.data[0] = cos(angle);
+ qexp.data[1] = omega(1) * sina_a;
+ qexp.data[2] = omega(2) * sina_a;
+ qexp.data[3] = omega(3) * sina_a;
+
+ return qexp;
+}
+
/** Write quaternion to a stream.
@param os Stream to write to.
@param q Quaternion to write.
FUNCTIONAL DESCRIPTION
--------------------------------------------------------------------------------
This class encapsulates the calculation of the derivatives of the state vectors
-UVW and PQR - the translational and rotational rates relative to the planet
+UVW and PQR - the translational and rotational rates relative to the planet
fixed frame. The derivatives relative to the inertial frame are also calculated
as a side effect. Also, the derivative of the attitude quaterion is also calculated.
[2] Richard E. McFarland, "A Standard Kinematic Model for Flight Simulation at
NASA-Ames", NASA CR-2497, January 1975
[3] Erin Catto, "Iterative Dynamics with Temporal Coherence", February 22, 2005
+[4] Mark Harris and Robert Lyle, "Spacecraft Gravitational Torques",
+ NASA SP-8024, May 1969
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
INCLUDES
namespace JSBSim {
-static const char *IdSrc = "$Id: FGAccelerations.cpp,v 1.8 2011/08/30 20:49:04 bcoconni Exp $";
+static const char *IdSrc = "$Id: FGAccelerations.cpp,v 1.13 2012/02/18 19:11:37 bcoconni Exp $";
static const char *IdHdr = ID_ACCELERATIONS;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Debug(0);
Name = "FGAccelerations";
gravType = gtWGS84;
-
+ gravTorque = false;
+
vPQRidot.InitMatrix();
vUVWidot.InitMatrix();
vGravAccel.InitMatrix();
if (FGModel::Run(Holding)) return true; // Fast return if we have nothing to do ...
if (Holding) return false;
- RunPreFunctions();
-
CalculatePQRdot(); // Angular rate derivative
CalculateUVWdot(); // Translational rate derivative
CalculateQuatdot(); // Angular orientation derivative
ResolveFrictionForces(in.DeltaT * rate); // Update rate derivatives with friction forces
- RunPostFunctions();
-
Debug(2);
return false;
}
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
// Compute body frame rotational accelerations based on the current body moments
//
-// vPQRdot is the derivative of the absolute angular velocity of the vehicle
-// (body rate with respect to the inertial frame), expressed in the body frame,
+// vPQRdot is the derivative of the absolute angular velocity of the vehicle
+// (body rate with respect to the ECEF frame), expressed in the body frame,
// where the derivative is taken in the body frame.
// J is the inertia matrix
// Jinv is the inverse inertia matrix
// vMoments is the moment vector in the body frame
// in.vPQRi is the total inertial angular velocity of the vehicle
// expressed in the body frame.
-// Reference: See Stevens and Lewis, "Aircraft Control and Simulation",
+// Reference: See Stevens and Lewis, "Aircraft Control and Simulation",
// Second edition (2004), eqn 1.5-16e (page 50)
void FGAccelerations::CalculatePQRdot(void)
{
+ if (gravTorque) {
+ // Compute the gravitational torque
+ // Reference: See Harris and Lyle "Spacecraft Gravitational Torques",
+ // NASA SP-8024 (1969) eqn (2) (page 7)
+ FGColumnVector3 R = in.Ti2b * in.vInertialPosition;
+ double invRadius = 1.0 / R.Magnitude();
+ R *= invRadius;
+ in.Moment += (3.0 * in.GAccel * invRadius) * (R * (in.J * R));
+ }
+
// Compute body frame rotational accelerations based on the current body
// moments and the total inertial angular velocity expressed in the body
// frame.
// Compute the quaternion orientation derivative
//
// vQtrndot is the quaternion derivative.
-// Reference: See Stevens and Lewis, "Aircraft Control and Simulation",
+// Reference: See Stevens and Lewis, "Aircraft Control and Simulation",
// Second edition (2004), eqn 1.5-16b (page 50)
void FGAccelerations::CalculateQuatdot(void)
// This set of calculations results in the body and inertial frame accelerations
// being computed.
// Compute body and inertial frames accelerations based on the current body
-// forces including centripetal and coriolis accelerations for the former.
+// forces including centripetal and Coriolis accelerations for the former.
// in.vOmegaPlanet is the Earth angular rate - expressed in the inertial frame -
// so it has to be transformed to the body frame. More completely,
// in.vOmegaPlanet is the rate of the ECEF frame relative to the Inertial
// in the body frame.
// in.vUVW is the vehicle velocity relative to the ECEF frame, expressed
// in the body frame.
-// Reference: See Stevens and Lewis, "Aircraft Control and Simulation",
+// Reference: See Stevens and Lewis, "Aircraft Control and Simulation",
// Second edition (2004), eqns 1.5-13 (pg 48) and 1.5-16d (page 50)
void FGAccelerations::CalculateUVWdot(void)
// Include Gravitation accel
switch (gravType) {
case gtStandard:
- vGravAccel = in.Tl2b * FGColumnVector3( 0.0, 0.0, in.GAccel );
+ {
+ double radius = in.vInertialPosition.Magnitude();
+ vGravAccel = -(in.GAccel / radius) * in.vInertialPosition;
+ }
break;
case gtWGS84:
- vGravAccel = in.Tec2b * in.J2Grav;
+ vGravAccel = in.Tec2i * in.J2Grav;
break;
}
- vUVWdot += vGravAccel;
- vUVWidot = in.Tb2i * (vBodyAccel + vGravAccel);
+ vUVWdot += in.Ti2b * vGravAccel;
+ vUVWidot = in.Tb2i * vBodyAccel + vGravAccel;
}
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
// Resolves the contact forces just before integrating the EOM.
// This routine is using Lagrange multipliers and the projected Gauss-Seidel
// (PGS) method.
-// Reference: See Erin Catto, "Iterative Dynamics with Temporal Coherence",
+// Reference: See Erin Catto, "Iterative Dynamics with Temporal Coherence",
// February 22, 2005
// In JSBSim there is only one rigid body (the aircraft) and there can be
// multiple points of contact between the aircraft and the ground. As a
-// consequence our matrix J*M^-1*J^T is not sparse and the algorithm described
-// in Catto's paper has been adapted accordingly.
+// consequence our matrix Jac*M^-1*Jac^T is not sparse and the algorithm
+// described in Catto's paper has been adapted accordingly.
// The friction forces are resolved in the body frame relative to the origin
// (Earth center).
// If no gears are in contact with the ground then return
if (!n) return;
- vector<double> a(n*n); // Will contain J*M^-1*J^T
+ vector<double> a(n*n); // Will contain Jac*M^-1*Jac^T
vector<double> rhs(n);
// Assemble the linear system of equations
FGColumnVector3 v2 = Jinv * multipliers[i]->MomentJacobian; // Should be J^-T but J is symmetric and so is J^-1
for (int j=0; j < i; j++)
- a[i*n+j] = a[j*n+i]; // Takes advantage of the symmetry of J^T*M^-1*J
+ a[i*n+j] = a[j*n+i]; // Takes advantage of the symmetry of Jac^T*M^-1*Jac
for (int j=i; j < n; j++)
a[i*n+j] = DotProduct(v1, multipliers[j]->ForceJacobian)
+ DotProduct(v2, multipliers[j]->MomentJacobian);
// Translation
vdot = vUVWdot;
- if (dt > 0.) // Zeroes out the relative movement between aircraft and ground
+ if (dt > 0.) // Zeroes out the relative movement between the aircraft and the ground
vdot += (in.vUVW - in.Tec2b * in.TerrainVelocity) / dt;
// Rotation
wdot = vPQRdot;
- if (dt > 0.) // Zeroes out the relative movement between aircraft and ground
+ if (dt > 0.) // Zeroes out the relative movement between the aircraft and the ground
wdot += (in.vPQR - in.Tec2b * in.TerrainAngularVel) / dt;
// Prepare the linear system for the Gauss-Seidel algorithm :
for (int i=0; i < n; i++) {
double lambda0 = multipliers[i]->value;
double dlambda = rhs[i];
-
+
for (int j=0; j < n; j++)
dlambda -= a[i*n+j]*multipliers[j]->value;
vPQRdot += omegadot;
vPQRidot += omegadot;
}
+
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
void FGAccelerations::InitializeDerivatives(void)
PropertyManager->Tie("accelerations/wdot-ft_sec2", this, eW, (PMF)&FGAccelerations::GetUVWdot);
PropertyManager->Tie("simulation/gravity-model", &gravType);
+ PropertyManager->Tie("simulation/gravitational-torque", &gravTorque);
PropertyManager->Tie("forces/fbx-total-lbs", this, eX, (PMF)&FGAccelerations::GetForces);
PropertyManager->Tie("forces/fby-total-lbs", this, eY, (PMF)&FGAccelerations::GetForces);
DEFINITIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#define ID_ACCELERATIONS "$Id: FGAccelerations.h,v 1.8 2011/10/31 14:54:41 bcoconni Exp $"
+#define ID_ACCELERATIONS "$Id: FGAccelerations.h,v 1.12 2012/02/19 14:07:27 bcoconni Exp $"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FORWARD DECLARATIONS
/** Handles the calculation of accelerations.
- -Calculate the angular accelerations
- -Calculate the translational accelerations
- -Calculate the angular rate
- -Calculate the translational velocity
+ - Calculate the angular accelerations
+ - Calculate the translational accelerations
+ - Calculate the angular rate
+ - Calculate the translational velocity
+
+ This class is collecting all the forces and the moments acting on the body
+ to calculate the corresponding accelerations according to Newton's second
+ law. This is also where the friction forces related to the ground reactions
+ are evaluated.
+
+ JSBSim provides several ways to calculate the influence of the gravity on
+ the vehicle. The different options can be selected via the following
+ properties :
+ @property simulation/gravity-model (read/write) Selects the gravity model.
+ Two options are available : 0 (Standard gravity assuming the Earth
+ is spherical) or 1 (WGS84 gravity taking the Earth oblateness into
+ account). WGS84 gravity is the default.
+ @property simulation/gravitational-torque (read/write) Enables/disables the
+ calculations of the gravitational torque on the vehicle. This is
+ mainly relevant for spacecrafts that are orbiting at low altitudes.
+ Gravitational torque calculations are disabled by default.
+
+ Special care is taken in the calculations to obtain maximum fidelity in
+ JSBSim results. In FGAccelerations, this is obtained by avoiding as much as
+ possible the transformations from one frame to another. As a consequence,
+ the frames in which the accelerations are primarily evaluated are dictated
+ by the frames in which FGPropagate resolves the equations of movement (the
+ ECI frame for the translations and the body frame for the rotations).
+
+ @see Mark Harris and Robert Lyle, "Spacecraft Gravitational Torques",
+ NASA SP-8024, May 1969
@author Jon S. Berndt, Mathias Froehlich, Bertrand Coconnier
- @version $Id: FGAccelerations.h,v 1.8 2011/10/31 14:54:41 bcoconni Exp $
+ @version $Id: FGAccelerations.h,v 1.12 2012/02/19 14:07:27 bcoconni Exp $
*/
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
/// Destructor
~FGAccelerations();
-
+
/// These define the indices use to select the gravitation models.
- enum eGravType {gtStandard, gtWGS84};
+ enum eGravType {
+ /// Evaluate gravity using Newton's classical formula assuming the Earth is spherical
+ gtStandard,
+ /// Evaluate gravity using WGS84 formulas that take the Earth oblateness into account
+ gtWGS84
+ };
/** Initializes the FGAccelerations class after instantiation and prior to first execution.
- The base class FGModel::InitModel is called first, initializing pointers to the
+ The base class FGModel::InitModel is called first, initializing pointers to the
other FGModel objects (and others). */
bool InitModel(void);
/** Runs the state propagation model; called by the Executive
Can pass in a value indicating if the executive is directing the simulation to Hold.
- @param Holding if true, the executive has been directed to hold the sim from
+ @param Holding if true, the executive has been directed to hold the sim from
advancing time. Some models may ignore this flag, such as the Input
model, which may need to be active to listen on a socket for the
"Resume" command to be given.
@return false if no error */
bool Run(bool Holding);
+ /** Retrieves the time derivative of the body orientation quaternion.
+ Retrieves the time derivative of the body orientation quaternion based on
+ the rate of change of the orientation between the body and the ECI frame.
+ The quaternion returned is represented by an FGQuaternion reference. The
+ quaternion is 1-based, so that the first element can be retrieved using
+ the "()" operator.
+ units rad/sec^2
+ @return The time derivative of the body orientation quaternion.
+ */
const FGQuaternion& GetQuaterniondot(void) const {return vQtrndot;}
/** Retrieves the body axis acceleration.
Retrieves the computed body axis accelerations based on the
applied forces and accounting for a rotating body frame.
- The vector returned is represented by an FGColumnVector reference. The vector
+ The vector returned is represented by an FGColumnVector3 reference. The vector
for the acceleration in Body frame is organized (Ax, Ay, Az). The vector
is 1-based, so that the first element can be retrieved using the "()" operator.
In other words, vUVWdot(1) is Ax. Various convenience enumerators are defined
*/
const FGColumnVector3& GetUVWdot(void) const { return vUVWdot; }
+ /** Retrieves the body axis acceleration in the ECI frame.
+ Retrieves the computed body axis accelerations based on the applied forces.
+ The ECI frame being an inertial frame this vector does not contain the
+ Coriolis and centripetal accelerations. The vector is expressed in the
+ Body frame.
+ The vector returned is represented by an FGColumnVector3 reference. The
+ vector for the acceleration in Body frame is organized (Aix, Aiy, Aiz). The
+ vector is 1-based, so that the first element can be retrieved using the
+ "()" operator. In other words, vUVWidot(1) is Aix. Various convenience
+ enumerators are defined in FGJSBBase. The relevant enumerators for the
+ vector returned by this call are, eX=1, eY=2, eZ=3.
+ units ft/sec^2
+ @return Body axis translational acceleration in ft/sec^2.
+ */
const FGColumnVector3& GetUVWidot(void) const { return vUVWidot; }
/** Retrieves the body axis angular acceleration vector.
Retrieves the body axis angular acceleration vector in rad/sec^2. The
- angular acceleration vector is determined from the applied forces and
+ angular acceleration vector is determined from the applied moments and
accounts for a rotating frame.
- The vector returned is represented by an FGColumnVector reference. The vector
+ The vector returned is represented by an FGColumnVector3 reference. The vector
for the angular acceleration in Body frame is organized (Pdot, Qdot, Rdot). The vector
is 1-based, so that the first element can be retrieved using the "()" operator.
In other words, vPQRdot(1) is Pdot. Various convenience enumerators are defined
@return The angular acceleration vector.
*/
const FGColumnVector3& GetPQRdot(void) const {return vPQRdot;}
-
+
+ /** Retrieves the axis angular acceleration vector in the ECI frame.
+ Retrieves the body axis angular acceleration vector measured in the ECI
+ frame and expressed in the body frame. The angular acceleration vector is
+ determined from the applied moments.
+ The vector returned is represented by an FGColumnVector3 reference. The
+ vector for the angular acceleration in Body frame is organized (Pidot,
+ Qidot, Ridot). The vector is 1-based, so that the first element can be
+ retrieved using the "()" operator. In other words, vPQRidot(1) is Pidot.
+ Various convenience enumerators are defined in FGJSBBase. The relevant
+ enumerators for the vector returned by this call are, eP=1, eQ=2, eR=3.
+ units rad/sec^2
+ @return The angular acceleration vector.
+ */
const FGColumnVector3& GetPQRidot(void) const {return vPQRidot;}
/** Retrieves a body frame acceleration component.
Retrieves a body frame acceleration component. The acceleration returned
- is extracted from the vUVWdot vector (an FGColumnVector). The vector for
+ is extracted from the vUVWdot vector (an FGColumnVector3). The vector for
the acceleration in Body frame is organized (Ax, Ay, Az). The vector is
1-based. In other words, GetUVWdot(1) returns Ax. Various convenience
enumerators are defined in FGJSBBase. The relevant enumerators for the
*/
double GetUVWdot(int idx) const { return vUVWdot(idx); }
+ /** Retrieves the acceleration resulting from the applied forces.
+ Retrieves the ratio of the sum of all forces applied on the craft to its
+ mass. This does include the friction forces but not the gravity.
+ The vector returned is represented by an FGColumnVector3 reference. The
+ vector for the acceleration in Body frame is organized (Ax, Ay, Az). The
+ vector is 1-based, so that the first element can be retrieved using the
+ "()" operator. In other words, vBodyAccel(1) is Ax. Various convenience
+ enumerators are defined in FGJSBBase. The relevant enumerators for the
+ vector returned by this call are, eX=1, eY=2, eZ=3.
+ units ft/sec^2
+ @return The acceleration resulting from the applied forces.
+ */
const FGColumnVector3& GetBodyAccel(void) const { return vBodyAccel; }
+
+ /** Retrieves a component of the acceleration resulting from the applied forces.
+ Retrieves a component of the ratio between the sum of all forces applied
+ on the craft to its mass. The value returned is extracted from the vBodyAccel
+ vector (an FGColumnVector3). The vector for the acceleration in Body frame
+ is organized (Ax, Ay, Az). The vector is 1-based. In other words,
+ GetBodyAccel(1) returns Ax. Various convenience enumerators are defined
+ in FGJSBBase. The relevant enumerators for the vector returned by this
+ call are, eX=1, eY=2, eZ=3.
+ units ft/sec^2
+ @param idx the index of the acceleration component desired (1-based).
+ @return The component of the acceleration resulting from the applied forces.
+ */
double GetBodyAccel(int idx) const { return vBodyAccel(idx); }
/** Retrieves a body frame angular acceleration component.
Retrieves a body frame angular acceleration component. The angular
acceleration returned is extracted from the vPQRdot vector (an
- FGColumnVector). The vector for the angular acceleration in Body frame
+ FGColumnVector3). The vector for the angular acceleration in Body frame
is organized (Pdot, Qdot, Rdot). The vector is 1-based. In other words,
GetPQRdot(1) returns Pdot (roll acceleration). Various convenience
enumerators are defined in FGJSBBase. The relevant enumerators for the
*/
double GetPQRdot(int axis) const {return vPQRdot(axis);}
+ /** Retrieves a component of the total moments applied on the body.
+ Retrieves a component of the total moments applied on the body. This does
+ include the moments generated by friction forces and the gravitational
+ torque (if the property \e simulation/gravitational-torque is set to true).
+ The vector for the total moments in the body frame is organized (Mx, My
+ , Mz). The vector is 1-based. In other words, GetMoments(1) returns Mx.
+ Various convenience enumerators are defined in FGJSBBase. The relevant
+ enumerators for the moments returned by this call are, eX=1, eY=2, eZ=3.
+ units lbs*ft
+ @param idx the index of the moments component desired (1-based).
+ @return The total moments applied on the body.
+ */
double GetMoments(int idx) const { return in.Moment(idx) + vFrictionMoments(idx); }
+
+ /** Retrieves the total forces applied on the body.
+ Retrieves the total forces applied on the body. This does include the
+ friction forces but not the gravity.
+ The vector for the total forces in the body frame is organized (Fx, Fy
+ , Fz). The vector is 1-based. In other words, GetForces(1) returns Fx.
+ Various convenience enumerators are defined in FGJSBBase. The relevant
+ enumerators for the forces returned by this call are, eX=1, eY=2, eZ=3.
+ units lbs
+ @param idx the index of the forces component desired (1-based).
+ @return The total forces applied on the body.
+ */
double GetForces(int idx) const { return in.Force(idx) + vFrictionForces(idx); }
+
+ /** Retrieves the ground moments applied on the body.
+ Retrieves the ground moments applied on the body. This does include the
+ ground normal reaction and friction moments.
+ The vector for the ground moments in the body frame is organized (Mx, My
+ , Mz). The vector is 1-based. In other words, GetGroundMoments(1) returns
+ Mx. Various convenience enumerators are defined in FGJSBBase. The relevant
+ enumerators for the moments returned by this call are, eX=1, eY=2, eZ=3.
+ units lbs*ft
+ @param idx the index of the moments component desired (1-based).
+ @return The ground moments applied on the body.
+ */
double GetGroundMoments(int idx) const { return in.GroundMoment(idx) + vFrictionMoments(idx); }
+
+ /** Retrieves the ground forces applied on the body.
+ Retrieves the ground forces applied on the body. This does include the
+ ground normal reaction and friction forces.
+ The vector for the total moments in the body frame is organized (Fx, Fy
+ , Fz). The vector is 1-based. In other words, GetGroundForces(1) returns
+ Fx. Various convenience enumerators are defined in FGJSBBase. The relevant
+ enumerators for the forces returned by this call are, eX=1, eY=2, eZ=3.
+ units lbs
+ @param idx the index of the forces component desired (1-based).
+ @return The ground forces applied on the body.
+ */
double GetGroundForces(int idx) const { return in.GroundForce(idx) + vFrictionForces(idx); }
+ /** Initializes the FGAccelerations class prior to a new execution.
+ Initializes the class prior to a new execution when the input data stored
+ in the Inputs structure have been set to their initial values.
+ */
void InitializeDerivatives(void);
- void DumpState(void);
-
struct Inputs {
+ /// The body inertia matrix expressed in the body frame
FGMatrix33 J;
+ /// The inverse of the inertia matrix J
FGMatrix33 Jinv;
+ /// Transformation matrix from the ECI to the Body frame
FGMatrix33 Ti2b;
+ /// Transformation matrix from the Body to the ECI frame
FGMatrix33 Tb2i;
+ /// Transformation matrix from the ECEF to the Body frame
FGMatrix33 Tec2b;
- FGMatrix33 Tl2b;
+ /// Transformation matrix from the ECEF to the ECI frame
+ FGMatrix33 Tec2i;
+ /// Orientation quaternion of the body with respect to the ECI frame
FGQuaternion qAttitudeECI;
+ /// Total moments applied to the body except friction and gravity (expressed in the body frame)
FGColumnVector3 Moment;
+ /// Moments generated by the ground normal reactions expressed in the body frame. Does not account for friction.
FGColumnVector3 GroundMoment;
+ /// Total forces applied to the body except friction and gravity (expressed in the body frame)
FGColumnVector3 Force;
+ /// Forces generated by the ground normal reactions expressed in the body frame. Does not account for friction.
FGColumnVector3 GroundForce;
+ /// Gravity intensity vector using WGS84 formulas (expressed in the ECEF frame).
FGColumnVector3 J2Grav;
+ /// Angular velocities of the body with respect to the ECI frame (expressed in the body frame).
FGColumnVector3 vPQRi;
+ /// Angular velocities of the body with respect to the local frame (expressed in the body frame).
FGColumnVector3 vPQR;
+ /// Velocities of the body with respect to the local frame (expressed in the body frame).
FGColumnVector3 vUVW;
+ /// Body position (X,Y,Z) measured in the ECI frame.
FGColumnVector3 vInertialPosition;
+ /// Earth rotating vector (expressed in the ECI frame).
FGColumnVector3 vOmegaPlanet;
+ /// Terrain velocities with respect to the local frame (expressed in the ECEF frame).
FGColumnVector3 TerrainVelocity;
+ /// Terrain angular velocities with respect to the local frame (expressed in the ECEF frame).
FGColumnVector3 TerrainAngularVel;
+ /// Time step
double DeltaT;
+ /// Body mass
double Mass;
+ /// Gravity intensity assuming the Earth is spherical
double GAccel;
+ /// List of Lagrange multipliers set by FGLGear for friction forces calculations.
std::vector<LagrangeMultiplier*> *MultipliersList;
} in;
FGColumnVector3 vFrictionMoments;
int gravType;
+ bool gravTorque;
void CalculatePQRdot(void);
void CalculateQuatdot(void);
namespace JSBSim {
-static const char *IdSrc = "$Id: FGAerodynamics.cpp,v 1.44 2011/10/23 14:23:13 jentron Exp $";
+static const char *IdSrc = "$Id: FGAerodynamics.cpp,v 1.45 2012/04/13 13:25:52 jberndt Exp $";
static const char *IdHdr = ID_AERODYNAMICS;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
AeroFunctionStrings += AeroFunctions[axis][sd]->GetName();
}
}
+
+ string FunctionStrings = FGModelFunctions::GetFunctionStrings(delimeter);
+
+ if (FunctionStrings.size() > 0) {
+ if (AeroFunctionStrings.size() > 0) {
+ AeroFunctionStrings += delimeter + FunctionStrings;
+ } else {
+ AeroFunctionStrings = FunctionStrings;
+ }
+ }
+
return AeroFunctionStrings;
}
}
}
+ string FunctionValues = FGModelFunctions::GetFunctionValues(delimeter);
+
+ if (FunctionValues.size() > 0) {
+ if (buf.str().size() > 0) {
+ buf << delimeter << FunctionValues;
+ } else {
+ buf << FunctionValues;
+ }
+ }
+
return buf.str();
}
namespace JSBSim {
-static const char *IdSrc = "$Id: FGAtmosphere.cpp,v 1.49 2011/09/11 11:36:04 bcoconni Exp $";
+static const char *IdSrc = "$Id: FGAtmosphere.cpp,v 1.51 2012/04/13 13:18:28 jberndt Exp $";
static const char *IdHdr = ID_ATMOSPHERE;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if (FGModel::Run(Holding)) return true;
if (Holding) return false;
- RunPreFunctions();
-
Calculate(in.altitudeASL);
- RunPostFunctions();
-
Debug(2);
return false;
}
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-void FGAtmosphere::SetPressureSL(double pressure, ePressure unit)
+void FGAtmosphere::SetPressureSL(ePressure unit, double pressure)
{
double press = ConvertToPSF(pressure, unit);
PropertyManager->Tie("atmosphere/a-fps", this, &FGAtmosphere::GetSoundSpeed);
PropertyManager->Tie("atmosphere/T-sl-R", this, &FGAtmosphere::GetTemperatureSL);
PropertyManager->Tie("atmosphere/rho-sl-slugs_ft3", this, &FGAtmosphere::GetDensitySL);
- PropertyManager->Tie("atmosphere/P-sl-psf", this, &FGAtmosphere::GetPressureSL);
+// PropertyManager->Tie("atmosphere/P-sl-psf", this, ePSF,
+// (PMFi)&FGAtmosphere::GetPressureSL,
+// (PMF)&FGAtmosphere::SetPressureSL);
PropertyManager->Tie("atmosphere/a-sl-fps", this, &FGAtmosphere::GetSoundSpeedSL);
PropertyManager->Tie("atmosphere/theta", this, &FGAtmosphere::GetTemperatureRatio);
PropertyManager->Tie("atmosphere/sigma", this, &FGAtmosphere::GetDensityRatio);
}
if (debug_lvl & 16) { // Sanity checking
}
- if (debug_lvl & 128) { //
+ if (debug_lvl & 128) { //
}
if (debug_lvl & 64) {
if (from == 0) { // Constructor
DEFINITIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#define ID_ATMOSPHERE "$Id: FGAtmosphere.h,v 1.29 2011/07/10 20:18:14 jberndt Exp $"
+#define ID_ATMOSPHERE "$Id: FGAtmosphere.h,v 1.30 2012/04/13 13:18:28 jberndt Exp $"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FORWARD DECLARATIONS
@property atmosphere/a-ratio
@author Jon Berndt
- @version $Id: FGAtmosphere.h,v 1.29 2011/07/10 20:18:14 jberndt Exp $
+ @version $Id: FGAtmosphere.h,v 1.30 2012/04/13 13:18:28 jberndt Exp $
*/
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
virtual double GetPressureRatio(void) const { return Pressure*rSLpressure; }
/** Sets the sea level pressure for modeling.
- @param pressure The pressure in the units specified (PSF by default).
+ @param pressure The pressure in the units specified.
@param unit the unit of measure that the specified pressure is
supplied in.*/
- virtual void SetPressureSL(double pressure, ePressure unit=ePSF);
+ virtual void SetPressureSL(ePressure unit, double pressure);
//@}
// *************************************************************************
namespace JSBSim {
-static const char *IdSrc = "$Id: FGAuxiliary.cpp,v 1.55 2011/11/12 18:59:11 bcoconni Exp $";
+static const char *IdSrc = "$Id: FGAuxiliary.cpp,v 1.56 2011/12/11 17:03:05 bcoconni Exp $";
static const char *IdHdr = ID_AUXILIARY;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
return true;
}
-
+
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FGAuxiliary::~FGAuxiliary()
if (FGModel::Run(Holding)) return true; // return true if error returned from base class
if (Holding) return false;
- RunPreFunctions();
-
-// Rotation
+ // Rotation
vEulerRates(eTht) = in.vPQR(eQ)*in.CosPhi - in.vPQR(eR)*in.SinPhi;
if (in.CosTht != 0.0) {
vEulerRates(ePhi) = in.vPQR(eP) + vEulerRates(ePsi)*in.SinTht;
}
-// Combine the wind speed with aircraft speed to obtain wind relative speed
+ // Combine the wind speed with aircraft speed to obtain wind relative speed
vAeroPQR = in.vPQR - in.TurbPQR;
vAeroUVW = in.vUVW - in.Tl2b * in.TotalWindNED;
vMachUVW(eW) = vAeroUVW(eW) / in.SoundSpeed;
double MachU2 = MachU * MachU;
-// Position
+ // Position
Vground = sqrt( in.vVel(eNorth)*in.vVel(eNorth) + in.vVel(eEast)*in.vVel(eEast) );
// When all models are executed calculate the distance from the initial point.
CalculateRelativePosition();
- RunPostFunctions();
-
return false;
}
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
void FGAuxiliary::CalculateRelativePosition(void) //ToDo: This belongs elsewhere - perhaps in FGPropagate or Exec
-{
+{
const double earth_radius_mt = in.ReferenceRadius*fttom;
lat_relative_position=(in.Latitude - FDMExec->GetIC()->GetLatitudeDegIC() *degtorad)*earth_radius_mt;
lon_relative_position=(in.Longitude - FDMExec->GetIC()->GetLongitudeDegIC()*degtorad)*earth_radius_mt;
namespace JSBSim {
-static const char *IdSrc = "$Id: FGGroundReactions.cpp,v 1.36 2011/08/21 15:13:22 bcoconni Exp $";
+static const char *IdSrc = "$Id: FGGroundReactions.cpp,v 1.39 2012/04/01 17:05:51 bcoconni Exp $";
static const char *IdHdr = ID_GROUNDREACTIONS;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FGModel::Load(el); // Perform base class Load
- in.vWhlBodyVec.resize(lGear.size());
-
for (unsigned int i=0; i<lGear.size();i++) lGear[i]->bind();
PostLoad(el, PropertyManager);
namespace JSBSim {
-static const char *IdSrc = "$Id: FGInertial.cpp,v 1.25 2011/10/31 14:54:41 bcoconni Exp $";
+static const char *IdSrc = "$Id: FGInertial.cpp,v 1.26 2011/12/11 17:03:05 bcoconni Exp $";
static const char *IdHdr = ID_INERTIAL;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
RadiusReference = 20925650.00; // Equatorial radius (WGS84)
C2_0 = -4.84165371736E-04; // WGS84 value for the C2,0 coefficient
J2 = 1.0826266836E-03; // WGS84 value for J2
- a = 20925646.3255; // WGS84 semimajor axis length in feet
+ a = 20925646.3255; // WGS84 semimajor axis length in feet
b = 20855486.5951; // WGS84 semiminor axis length in feet
// Lunar defaults
RadiusReference = 5702559.05; // Equatorial radius
C2_0 = 0; // value for the C2,0 coefficient
J2 = 2.033542482111609E-4; // value for J2
- a = 5702559.05; // semimajor axis length in feet
+ a = 5702559.05; // semimajor axis length in feet
b = 5695439.63; // semiminor axis length in feet
*/
if (FGModel::Run(Holding)) return true;
if (Holding) return false;
- RunPreFunctions();
-
// Gravitation accel
gAccel = GetGAccel(in.Radius);
- RunPostFunctions();
-
return false;
}
namespace JSBSim {
-static const char *IdSrc = "$Id: FGInput.cpp,v 1.21 2011/05/20 03:18:36 jberndt Exp $";
+static const char *IdSrc = "$Id: FGInput.cpp,v 1.22 2012/01/21 16:46:09 jberndt Exp $";
static const char *IdHdr = ID_INPUT;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FDMExec->Resume();
socket->Reply("");
+
+ } else if (command == "iterate") { // ITERATE
+
+ int argumentInt;
+ istringstream (argument) >> argumentInt;
+ if (argument.size() == 0) {
+ socket->Reply("No argument supplied for number of iterations.\n");
+ break;
+ }
+ if ( !(argumentInt > 0) ){
+ socket->Reply("Required argument must be a positive Integer.\n");
+ break;
+ }
+ FDMExec->EnableIncrementThenHold( argumentInt );
+ FDMExec->Resume();
+ socket->Reply("");
} else if (command == "quit") { // QUIT
GLOBAL DATA
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-static const char *IdSrc = "$Id: FGLGear.cpp,v 1.92 2011/11/10 12:06:14 jberndt Exp $";
+static const char *IdSrc = "$Id: FGLGear.cpp,v 1.100 2012/04/01 17:05:51 bcoconni Exp $";
static const char *IdHdr = ID_LGEAR;
// Body To Structural (body frame is rotated 180 deg about Y and lengths are given in
// ft instead of inches)
const FGMatrix33 FGLGear::Tb2s(-1./inchtoft, 0., 0., 0., 1./inchtoft, 0., 0., 0., -1./inchtoft);
+const FGMatrix33 FGLGear::Ts2b(-inchtoft, 0., 0., 0., inchtoft, 0., 0., 0., -inchtoft);
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
CLASS IMPLEMENTATION
Castered(false),
StaticFriction(false)
{
- Element *force_table=0;
- Element *dampCoeff=0;
- Element *dampCoeffRebound=0;
- string force_type="";
-
kSpring = bDamp = bDampRebound = dynamicFCoeff = staticFCoeff = rollingFCoeff = maxSteerAngle = 0;
- sSteerType = sBrakeGroup = sSteerType = "";
- isRetractable = 0;
+ isRetractable = false;
eDampType = dtLinear;
eDampTypeRebound = dtLinear;
name = el->GetAttributeValue("name");
- sContactType = el->GetAttributeValue("type");
+ string sContactType = el->GetAttributeValue("type");
if (sContactType == "BOGEY") {
eContactType = ctBOGEY;
} else if (sContactType == "STRUCTURE") {
eContactType = ctSTRUCTURE;
}
- // Default values for structural contact points
+ // Default values for structural contact points
if (eContactType == ctSTRUCTURE) {
kSpring = in.EmptyWeight;
bDamp = kSpring;
if (el->FindElement("spring_coeff"))
kSpring = el->FindElementValueAsNumberConvertTo("spring_coeff", "LBS/FT");
if (el->FindElement("damping_coeff")) {
- dampCoeff = el->FindElement("damping_coeff");
+ Element* dampCoeff = el->FindElement("damping_coeff");
if (dampCoeff->GetAttributeValue("type") == "SQUARE") {
eDampType = dtSquare;
bDamp = el->FindElementValueAsNumberConvertTo("damping_coeff", "LBS/FT2/SEC2");
}
if (el->FindElement("damping_coeff_rebound")) {
- dampCoeffRebound = el->FindElement("damping_coeff_rebound");
+ Element* dampCoeffRebound = el->FindElement("damping_coeff_rebound");
if (dampCoeffRebound->GetAttributeValue("type") == "SQUARE") {
eDampTypeRebound = dtSquare;
bDampRebound = el->FindElementValueAsNumberConvertTo("damping_coeff_rebound", "LBS/FT2/SEC2");
staticFCoeff = el->FindElementValueAsNumber("static_friction");
if (el->FindElement("rolling_friction"))
rollingFCoeff = el->FindElementValueAsNumber("rolling_friction");
- if (el->FindElement("max_steer"))
- maxSteerAngle = el->FindElementValueAsNumberConvertTo("max_steer", "DEG");
if (el->FindElement("retractable"))
isRetractable = ((unsigned int)el->FindElementValueAsNumber("retractable"))>0.0?true:false;
+ if (el->FindElement("max_steer"))
+ maxSteerAngle = el->FindElementValueAsNumberConvertTo("max_steer", "DEG");
+
+ if (maxSteerAngle == 360) {
+ eSteerType = stCaster;
+ Castered = true;
+ }
+ else if (maxSteerAngle == 0.0) {
+ eSteerType = stFixed;
+ }
+ else
+ eSteerType = stSteer;
+
GroundReactions = fdmex->GetGroundReactions();
PropertyManager = fdmex->GetPropertyManager();
ForceY_Table = 0;
- force_table = el->FindElement("table");
+ Element* force_table = el->FindElement("table");
while (force_table) {
- force_type = force_table->GetAttributeValue("type");
+ string force_type = force_table->GetAttributeValue("type");
if (force_type == "CORNERING_COEFF") {
ForceY_Table = new FGTable(PropertyManager, force_table);
+ break;
} else {
cerr << "Undefined force table for " << name << " contact point" << endl;
}
force_table = el->FindNextElement("table");
}
- sBrakeGroup = el->FindElementValue("brake_group");
-
- if (maxSteerAngle == 360) sSteerType = "CASTERED";
- else if (maxSteerAngle == 0.0) sSteerType = "FIXED";
- else sSteerType = "STEERABLE";
-
Element* element = el->FindElement("location");
if (element) vXYZn = element->FindElementTripletConvertTo("IN");
else {cerr << "No location given for contact " << name << endl; exit(-1);}
element = el->FindElement("orientation");
if (element && (eContactType == ctBOGEY)) {
- vGearOrient = element->FindElementTripletConvertTo("RAD");
-
- double cp,sp,cr,sr,cy,sy;
+ FGQuaternion quatFromEuler(element->FindElementTripletConvertTo("RAD"));
- cp=cos(vGearOrient(ePitch)); sp=sin(vGearOrient(ePitch));
- cr=cos(vGearOrient(eRoll)); sr=sin(vGearOrient(eRoll));
- cy=cos(vGearOrient(eYaw)); sy=sin(vGearOrient(eYaw));
-
- mTGear(1,1) = cp*cy;
- mTGear(2,1) = cp*sy;
- mTGear(3,1) = -sp;
-
- mTGear(1,2) = sr*sp*cy - cr*sy;
- mTGear(2,2) = sr*sp*sy + cr*cy;
- mTGear(3,2) = sr*cp;
-
- mTGear(1,3) = cr*sp*cy + sr*sy;
- mTGear(2,3) = cr*sp*sy - sr*cy;
- mTGear(3,3) = cr*cp;
+ mTGear = quatFromEuler.GetT();
}
else {
mTGear(1,1) = 1.;
mTGear(3,3) = 1.;
}
+ string sBrakeGroup = el->FindElementValue("brake_group");
+
if (sBrakeGroup == "LEFT" ) eBrakeGrp = bgLeft;
else if (sBrakeGroup == "RIGHT" ) eBrakeGrp = bgRight;
else if (sBrakeGroup == "CENTER") eBrakeGrp = bgCenter;
- else if (sBrakeGroup == "NOSE" ) eBrakeGrp = bgNose;
- else if (sBrakeGroup == "TAIL" ) eBrakeGrp = bgTail;
+ else if (sBrakeGroup == "NOSE" ) eBrakeGrp = bgCenter; // Nose brake is not supported by FGFCS
+ else if (sBrakeGroup == "TAIL" ) eBrakeGrp = bgCenter; // Tail brake is not supported by FGFCS
else if (sBrakeGroup == "NONE" ) eBrakeGrp = bgNone;
- else if (sBrakeGroup.empty() ) {eBrakeGrp = bgNone;
- sBrakeGroup = "NONE (defaulted)";}
+ else if (sBrakeGroup.empty() ) eBrakeGrp = bgNone;
else {
cerr << "Improper braking group specification in config file: "
<< sBrakeGroup << " is undefined." << endl;
}
- if (sSteerType == "STEERABLE") eSteerType = stSteer;
- else if (sSteerType == "FIXED" ) eSteerType = stFixed;
- else if (sSteerType == "CASTERED" ) {eSteerType = stCaster; Castered = true;}
- else if (sSteerType.empty() ) {eSteerType = stFixed;
- sSteerType = "FIXED (defaulted)";}
- else {
- cerr << "Improper steering type specification in config file: "
- << sSteerType << " is undefined." << endl;
- }
-
- GearUp = false;
- GearDown = true;
GearPos = 1.0;
useFCSGearPos = false;
- Servicable = true;
// Add some AI here to determine if gear is located properly according to its
// brake group type ??
compressLength = 0.0;
compressSpeed = 0.0;
- brakePct = 0.0;
maxCompLen = 0.0;
WheelSlip = 0.0;
- TirePressureNorm = 1.0;
// Set Pacejka terms
const FGColumnVector3& FGLGear::GetBodyForces(void)
{
+ double gearPos = 1.0;
double t = fdmex->GetSimTime();
vFn.InitMatrix();
- if (isRetractable) ComputeRetractionState();
+ if (isRetractable) gearPos = GetGearUnitPos();
- if (GearDown) {
- FGColumnVector3 terrainVel, dummy;
-
- vLocalGear = in.Tb2l * in.vWhlBodyVec[GearNumber]; // Get local frame wheel location
+ if (gearPos > 0.99) { // Gear DOWN
+ FGColumnVector3 normal, terrainVel, dummy;
+ FGLocation gearLoc, contact;
+ FGColumnVector3 vWhlBodyVec = Ts2b * (vXYZn - in.vXYZcg);
+ vLocalGear = in.Tb2l * vWhlBodyVec; // Get local frame wheel location
gearLoc = in.Location.LocalToLocation(vLocalGear);
+
// Compute the height of the theoretical location of the wheel (if strut is
// not compressed) with respect to the ground level
double height = gearLoc.GetContactPoint(t, contact, normal, terrainVel, dummy);
- vGroundNormal = in.Tec2b * normal;
- // The height returned above is the AGL and is expressed in the Z direction
- // of the ECEF coordinate frame. We now need to transform this height in
- // actual compression of the strut (BOGEY) of in the normal direction to the
- // ground (STRUCTURE)
- double normalZ = (in.Tec2l*normal)(eZ);
- double LGearProj = -(mTGear.Transposed() * vGroundNormal)(eZ);
-
- switch (eContactType) {
- case ctBOGEY:
- compressLength = LGearProj > 0.0 ? height * normalZ / LGearProj : 0.0;
- break;
- case ctSTRUCTURE:
- compressLength = height * normalZ / DotProduct(normal, normal);
- break;
- }
-
- if (compressLength > 0.00) {
+ if (height < 0.0) {
WOW = true;
+ vGroundNormal = in.Tec2b * normal;
+
+ // The height returned by GetGroundCallback() is the AGL and is expressed
+ // in the Z direction of the local coordinate frame. We now need to transform
+ // this height in actual compression of the strut (BOGEY) or in the normal
+ // direction to the ground (STRUCTURE)
+ double normalZ = (in.Tec2l*normal)(eZ);
+ double LGearProj = -(mTGear.Transposed() * vGroundNormal)(eZ);
+ FGColumnVector3 vWhlDisplVec;
// The following equations use the vector to the tire contact patch
// including the strut compression.
- FGColumnVector3 vWhlDisplVec;
-
switch(eContactType) {
case ctBOGEY:
+ compressLength = LGearProj > 0.0 ? height * normalZ / LGearProj : 0.0;
vWhlDisplVec = mTGear * FGColumnVector3(0., 0., -compressLength);
break;
case ctSTRUCTURE:
+ compressLength = height * normalZ / DotProduct(normal, normal);
vWhlDisplVec = compressLength * vGroundNormal;
break;
}
- FGColumnVector3 vWhlContactVec = in.vWhlBodyVec[GearNumber] + vWhlDisplVec;
+ FGColumnVector3 vWhlContactVec = vWhlBodyVec + vWhlDisplVec;
vActingXYZn = vXYZn + Tb2s * vWhlDisplVec;
FGColumnVector3 vBodyWhlVel = in.PQR * vWhlContactVec;
vBodyWhlVel += in.UVW - in.Tec2b * terrainVel;
InitializeReporting();
ComputeSteeringAngle();
- ComputeGroundCoordSys();
+ ComputeGroundFrame();
- vLocalWhlVel = Transform().Transposed() * vBodyWhlVel;
+ vGroundWhlVel = mT.Transposed() * vBodyWhlVel;
if (fdmex->GetTrimStatus())
compressSpeed = 0.0; // Steady state is sought during trimming
else {
- compressSpeed = -vLocalWhlVel(eX);
- if (eContactType == ctBOGEY) compressSpeed /= LGearProj;
+ compressSpeed = -vGroundWhlVel(eZ);
+ if (eContactType == ctBOGEY)
+ compressSpeed /= LGearProj;
}
ComputeVerticalStrutForce();
WheelSlip = 0.0;
StrutForce = 0.0;
+ LMultiplier[ftRoll].value = 0.0;
+ LMultiplier[ftSide].value = 0.0;
+ LMultiplier[ftDynamic].value = 0.0;
+
// Let wheel spin down slowly
vWhlVelVec(eX) -= 13.0 * in.TotalDeltaT;
if (vWhlVelVec(eX) < 0.0) vWhlVelVec(eX) = 0.0;
// Return to neutral position between 1.0 and 0.8 gear pos.
- SteerAngle *= max(GetGearUnitPos()-0.8, 0.0)/0.2;
+ SteerAngle *= max(gearPos-0.8, 0.0)/0.2;
ResetReporting();
}
}
+ else if (gearPos < 0.01) { // Gear UP
+ WOW = false;
+ vWhlVelVec.InitMatrix();
+ }
if (!fdmex->GetTrimStatus()) {
ReportTakeoffOrLanding();
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
// Build a local "ground" coordinate system defined by
-// eX : normal to the ground
-// eY : projection of the rolling direction on the ground
-// eZ : projection of the sliping direction on the ground
+// eX : projection of the rolling direction on the ground
+// eY : projection of the sliping direction on the ground
+// eZ : normal to the ground
-void FGLGear::ComputeGroundCoordSys(void)
+void FGLGear::ComputeGroundFrame(void)
{
- // Euler angles are built up to create a local frame to describe the forces
- // applied to the gear by the ground. Here pitch, yaw and roll do not have
- // any physical meaning. It is just a convenient notation.
- // First, "pitch" and "yaw" are determined in order to align eX with the
- // ground normal.
- if (vGroundNormal(eZ) < -1.0)
- vOrient(ePitch) = 0.5*M_PI;
- else if (1.0 < vGroundNormal(eZ))
- vOrient(ePitch) = -0.5*M_PI;
- else
- vOrient(ePitch) = asin(-vGroundNormal(eZ));
-
- if (fabs(vOrient(ePitch)) == 0.5*M_PI)
- vOrient(eYaw) = 0.;
- else
- vOrient(eYaw) = atan2(vGroundNormal(eY), vGroundNormal(eX));
-
- vOrient(eRoll) = 0.;
- UpdateCustomTransformMatrix();
-
- if (eContactType == ctBOGEY) {
- // In the case of a bogey, the third angle "roll" is used to align the axis eY and eZ
- // to the rolling and sliping direction respectively.
- FGColumnVector3 updatedRollingAxis = Transform().Transposed() * mTGear
- * FGColumnVector3(-sin(SteerAngle), cos(SteerAngle), 0.);
-
- vOrient(eRoll) = atan2(updatedRollingAxis(eY), -updatedRollingAxis(eZ));
- UpdateCustomTransformMatrix();
- }
-}
-
-//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-
-void FGLGear::ComputeRetractionState(void)
-{
- double gearPos = GetGearUnitPos();
- if (gearPos < 0.01) {
- GearUp = true;
- WOW = false;
- GearDown = false;
- vWhlVelVec.InitMatrix();
- } else if (gearPos > 0.99) {
- GearDown = true;
- GearUp = false;
- } else {
- GearUp = false;
- GearDown = false;
- }
+ FGColumnVector3 roll = mTGear * FGColumnVector3(cos(SteerAngle), sin(SteerAngle), 0.);
+ FGColumnVector3 side = vGroundNormal * roll;
+
+ roll -= DotProduct(roll, vGroundNormal) * vGroundNormal;
+ roll.Normalize();
+ side.Normalize();
+
+ mT(eX,eX) = roll(eX);
+ mT(eY,eX) = roll(eY);
+ mT(eZ,eX) = roll(eZ);
+ mT(eX,eY) = side(eX);
+ mT(eY,eY) = side(eY);
+ mT(eZ,eY) = side(eZ);
+ mT(eX,eZ) = vGroundNormal(eX);
+ mT(eY,eZ) = vGroundNormal(eY);
+ mT(eZ,eZ) = vGroundNormal(eZ);
}
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
void FGLGear::ComputeSlipAngle(void)
{
// Check that the speed is non-null otherwise use the current angle
- if (vLocalWhlVel.Magnitude(eY,eZ) > 1E-3)
- WheelSlip = -atan2(vLocalWhlVel(eZ), fabs(vLocalWhlVel(eY)))*radtodeg;
+ if (vGroundWhlVel.Magnitude(eX,eY) > 1E-3)
+ WheelSlip = -atan2(vGroundWhlVel(eY), fabs(vGroundWhlVel(eX)))*radtodeg;
}
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
void FGLGear::ComputeBrakeForceCoefficient(void)
{
- switch (eBrakeGrp) {
- case bgLeft:
- BrakeFCoeff = ( rollingFCoeff * (1.0 - in.BrakePos[bgLeft]) +
- staticFCoeff * in.BrakePos[bgLeft] );
- break;
- case bgRight:
- BrakeFCoeff = ( rollingFCoeff * (1.0 - in.BrakePos[bgRight]) +
- staticFCoeff * in.BrakePos[bgRight] );
- break;
- case bgCenter:
- BrakeFCoeff = ( rollingFCoeff * (1.0 - in.BrakePos[bgCenter]) +
- staticFCoeff * in.BrakePos[bgCenter] );
- break;
- case bgNose:
- BrakeFCoeff = ( rollingFCoeff * (1.0 - in.BrakePos[bgCenter]) +
- staticFCoeff * in.BrakePos[bgCenter] );
- break;
- case bgTail:
- BrakeFCoeff = ( rollingFCoeff * (1.0 - in.BrakePos[bgCenter]) +
- staticFCoeff * in.BrakePos[bgCenter] );
- break;
- case bgNone:
- BrakeFCoeff = rollingFCoeff;
- break;
- default:
- cerr << "Improper brake group membership detected for this gear." << endl;
- break;
- }
+ BrakeFCoeff = rollingFCoeff;
+
+ if (eBrakeGrp != bgNone)
+ BrakeFCoeff += in.BrakePos[eBrakeGrp] * (staticFCoeff - rollingFCoeff);
}
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if (compressSpeed >= 0.0) {
- if (eDampType == dtLinear) dampForce = -compressSpeed * bDamp;
- else dampForce = -compressSpeed * compressSpeed * bDamp;
+ if (eDampType == dtLinear)
+ dampForce = -compressSpeed * bDamp;
+ else
+ dampForce = -compressSpeed * compressSpeed * bDamp;
} else {
switch (eContactType) {
case ctBOGEY:
// Project back the strut force in the local coordinate frame of the ground
- vFn(eX) = StrutForce / (mTGear.Transposed()*vGroundNormal)(eZ);
+ vFn(eZ) = StrutForce / (mTGear.Transposed()*vGroundNormal)(eZ);
break;
case ctSTRUCTURE:
- vFn(eX) = -StrutForce;
+ vFn(eZ) = -StrutForce;
break;
}
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-double FGLGear::GetGearUnitPos(void)
+double FGLGear::GetGearUnitPos(void) const
{
// hack to provide backward compatibility to gear/gear-pos-norm property
if( useFCSGearPos || in.FCSGearPos != 1.0 ) {
// When the point of contact is moving, dynamic friction is used
// This type of friction is limited to ctSTRUCTURE elements because their
// friction coefficient is the same in every directions
- if ((eContactType == ctSTRUCTURE) && (vLocalWhlVel.Magnitude(eY,eZ) > 1E-3)) {
- FGColumnVector3 velocityDirection = vLocalWhlVel;
+ if ((eContactType == ctSTRUCTURE) && (vGroundWhlVel.Magnitude(eX,eY) > 1E-3)) {
+
+ FGColumnVector3 velocityDirection = vGroundWhlVel;
StaticFriction = false;
- velocityDirection(eX) = 0.;
+ velocityDirection(eZ) = 0.;
velocityDirection.Normalize();
- LMultiplier[ftDynamic].ForceJacobian = Transform()*velocityDirection;
+ LMultiplier[ftDynamic].ForceJacobian = mT * velocityDirection;
LMultiplier[ftDynamic].MomentJacobian = vWhlContactVec * LMultiplier[ftDynamic].ForceJacobian;
LMultiplier[ftDynamic].Max = 0.;
- LMultiplier[ftDynamic].Min = -fabs(dynamicFCoeff * vFn(eX));
+ LMultiplier[ftDynamic].Min = -fabs(dynamicFCoeff * vFn(eZ));
// The Lagrange multiplier value obtained from the previous iteration is kept
// This is supposed to accelerate the convergence of the projected Gauss-Seidel
// This cannot be handled properly by the so-called "dynamic friction".
StaticFriction = true;
- LMultiplier[ftRoll].ForceJacobian = Transform()*FGColumnVector3(0.,1.,0.);
- LMultiplier[ftSide].ForceJacobian = Transform()*FGColumnVector3(0.,0.,1.);
+ LMultiplier[ftRoll].ForceJacobian = mT * FGColumnVector3(1.,0.,0.);
+ LMultiplier[ftSide].ForceJacobian = mT * FGColumnVector3(0.,1.,0.);
LMultiplier[ftRoll].MomentJacobian = vWhlContactVec * LMultiplier[ftRoll].ForceJacobian;
LMultiplier[ftSide].MomentJacobian = vWhlContactVec * LMultiplier[ftSide].ForceJacobian;
switch(eContactType) {
case ctBOGEY:
- LMultiplier[ftRoll].Max = fabs(BrakeFCoeff * vFn(eX));
- LMultiplier[ftSide].Max = fabs(FCoeff * vFn(eX));
+ LMultiplier[ftRoll].Max = fabs(BrakeFCoeff * vFn(eZ));
+ LMultiplier[ftSide].Max = fabs(FCoeff * vFn(eZ));
break;
case ctSTRUCTURE:
- LMultiplier[ftRoll].Max = fabs(staticFCoeff * vFn(eX));
- LMultiplier[ftSide].Max = fabs(staticFCoeff * vFn(eX));
+ LMultiplier[ftRoll].Max = fabs(staticFCoeff * vFn(eZ));
+ LMultiplier[ftSide].Max = LMultiplier[ftRoll].Max;
break;
}
void FGLGear::UpdateForces(void)
{
if (StaticFriction) {
- vFn(eY) = LMultiplier[ftRoll].value;
- vFn(eZ) = LMultiplier[ftSide].value;
+ vFn(eX) = LMultiplier[ftRoll].value;
+ vFn(eY) = LMultiplier[ftSide].value;
+ }
+ else {
+ FGColumnVector3 forceDir = mT.Transposed() * LMultiplier[ftDynamic].ForceJacobian;
+ vFn(eX) = LMultiplier[ftDynamic].value * forceDir(eX);
+ vFn(eY) = LMultiplier[ftDynamic].value * forceDir(eY);
}
- else
- vFn += LMultiplier[ftDynamic].value * (Transform ().Transposed() * LMultiplier[ftDynamic].ForceJacobian);
}
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
void FGLGear::Debug(int from)
{
+ static const char* sSteerType[] = {"STEERABLE", "FIXED", "CASTERED" };
+ static const char* sBrakeGroup[] = {"NONE", "LEFT", "RIGHT", "CENTER", "NOSE", "TAIL"};
+ static const char* sContactType[] = {"BOGEY", "STRUCTURE" };
+
if (debug_lvl <= 0) return;
if (debug_lvl & 1) { // Standard console startup message output
if (from == 0) { // Constructor - loading and initialization
- cout << " " << sContactType << " " << name << endl;
+ cout << " " << sContactType[eContactType] << " " << name << endl;
cout << " Location: " << vXYZn << endl;
cout << " Spring Constant: " << kSpring << endl;
if (eDampTypeRebound == dtLinear)
cout << " Rebound Damping Constant: " << bDampRebound << " (linear)" << endl;
- else
+ else
cout << " Rebound Damping Constant: " << bDampRebound << " (square law)" << endl;
cout << " Dynamic Friction: " << dynamicFCoeff << endl;
cout << " Static Friction: " << staticFCoeff << endl;
if (eContactType == ctBOGEY) {
cout << " Rolling Friction: " << rollingFCoeff << endl;
- cout << " Steering Type: " << sSteerType << endl;
- cout << " Grouping: " << sBrakeGroup << endl;
+ cout << " Steering Type: " << sSteerType[eSteerType] << endl;
+ cout << " Grouping: " << sBrakeGroup[eBrakeGrp] << endl;
cout << " Max Steer Angle: " << maxSteerAngle << endl;
cout << " Retractable: " << isRetractable << endl;
}
}
} // namespace JSBSim
-
#include "models/propulsion/FGForce.h"
#include "math/FGColumnVector3.h"
-#include "math/FGMatrix33.h"
#include "math/LagrangeMultiplier.h"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
DEFINITIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#define ID_LGEAR "$Id: FGLGear.h,v 1.48 2011/10/31 14:54:41 bcoconni Exp $"
+#define ID_LGEAR "$Id: FGLGear.h,v 1.54 2012/04/01 17:05:51 bcoconni Exp $"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FORWARD DECLARATIONS
</contact>
@endcode
@author Jon S. Berndt
- @version $Id: FGLGear.h,v 1.48 2011/10/31 14:54:41 bcoconni Exp $
+ @version $Id: FGLGear.h,v 1.54 2012/04/01 17:05:51 bcoconni Exp $
@see Richard E. McFarland, "A Standard Kinematic Model for Flight Simulation at
NASA-Ames", NASA CR-2497, January 1975
@see Barnes W. McCormick, "Aerodynamics, Aeronautics, and Flight Mechanics",
FGMatrix33 Tec2b;
FGColumnVector3 PQR;
FGColumnVector3 UVW;
+ FGColumnVector3 vXYZcg; // CG coordinates expressed in the structural frame
FGLocation Location;
std::vector <double> SteerPosDeg;
std::vector <double> BrakePos;
- std::vector <FGColumnVector3> vWhlBodyVec;
double FCSGearPos;
double EmptyWeight;
};
const FGColumnVector3& GetBodyForces(void);
/// Gets the location of the gear in Body axes
- const FGColumnVector3& GetBodyLocation(void) const { return in.vWhlBodyVec[GearNumber]; }
- double GetBodyLocation(int idx) const { return in.vWhlBodyVec[GearNumber](idx); }
+ FGColumnVector3 GetBodyLocation(void) const {
+ return Ts2b * (vXYZn - in.vXYZcg);
+ }
+ double GetBodyLocation(int idx) const {
+ FGColumnVector3 vWhlBodyVec = Ts2b * (vXYZn - in.vXYZcg);
+ return vWhlBodyVec(idx);
+ }
const FGColumnVector3& GetLocalGear(void) const { return vLocalGear; }
double GetLocalGear(int idx) const { return vLocalGear(idx); }
double GetCompVel(void) const {return compressSpeed; }
/// Gets the gear compression force in pounds
double GetCompForce(void) const {return StrutForce; }
- double GetBrakeFCoeff(void) const {return BrakeFCoeff;}
-
- /// Gets the current normalized tire pressure
- double GetTirePressure(void) const { return TirePressureNorm; }
- /// Sets the new normalized tire pressure
- void SetTirePressure(double p) { TirePressureNorm = p; }
-
- /// Sets the brake value in percent (0 - 100)
- void SetBrake(double bp) {brakePct = bp;}
/// Sets the weight-on-wheels flag.
void SetWOW(bool wow) {WOW = wow;}
bool GetSteerable(void) const { return eSteerType != stFixed; }
bool GetRetractable(void) const { return isRetractable; }
- bool GetGearUnitUp(void) const { return GearUp; }
- bool GetGearUnitDown(void) const { return GearDown; }
+ bool GetGearUnitUp(void) const { return isRetractable ? (GetGearUnitPos() < 0.01) : false; }
+ bool GetGearUnitDown(void) const { return isRetractable ? (GetGearUnitPos() > 0.99) : true; }
+
double GetWheelRollForce(void) {
UpdateForces();
FGColumnVector3 vForce = mTGear.Transposed() * FGForce::GetBodyForces();
double GetWheelSlipAngle(void) const { return WheelSlip; }
double GetWheelVel(int axis) const { return vWhlVelVec(axis);}
bool IsBogey(void) const { return (eContactType == ctBOGEY);}
- double GetGearUnitPos(void);
+ double GetGearUnitPos(void) const;
double GetSteerAngleDeg(void) const { return radtodeg*SteerAngle; }
const struct Inputs& in;
private:
int GearNumber;
- static const FGMatrix33 Tb2s;
+ static const FGMatrix33 Tb2s, Ts2b;
FGMatrix33 mTGear;
- FGColumnVector3 vGearOrient;
FGColumnVector3 vLocalGear;
- FGColumnVector3 vWhlVelVec, vLocalWhlVel; // Velocity of this wheel
- FGColumnVector3 normal, vGroundNormal;
- FGLocation contact, gearLoc;
+ FGColumnVector3 vWhlVelVec, vGroundWhlVel; // Velocity of this wheel
+ FGColumnVector3 vGroundNormal;
FGTable *ForceY_Table;
double SteerAngle;
double kSpring;
double compressSpeed;
double staticFCoeff, dynamicFCoeff, rollingFCoeff;
double Stiffness, Shape, Peak, Curvature; // Pacejka factors
- double brakePct;
double BrakeFCoeff;
double maxCompLen;
double SinkRate;
double MaximumStrutTravel;
double FCoeff;
double WheelSlip;
- double TirePressureNorm;
double GearPos;
- bool useFCSGearPos;
bool WOW;
bool lastWOW;
bool FirstContact;
bool TakeoffReported;
bool ReportEnable;
bool isRetractable;
- bool GearUp, GearDown;
- bool Servicable;
bool Castered;
bool StaticFriction;
std::string name;
- std::string sSteerType;
- std::string sBrakeGroup;
- std::string sRetractable;
- std::string sContactType;
BrakeGroup eBrakeGrp;
ContactType eContactType;
FGGroundReactions* GroundReactions;
FGPropertyManager* PropertyManager;
- void ComputeRetractionState(void);
+ mutable bool useFCSGearPos;
+
void ComputeBrakeForceCoefficient(void);
void ComputeSteeringAngle(void);
void ComputeSlipAngle(void);
void ComputeSideForceCoefficient(void);
void ComputeVerticalStrutForce(void);
- void ComputeGroundCoordSys(void);
+ void ComputeGroundFrame(void);
void ComputeJacobian(const FGColumnVector3& vWhlContactVec);
void UpdateForces(void);
void CrashDetect(void);
namespace JSBSim {
-static const char *IdSrc = "$Id: FGOutput.cpp,v 1.63 2011/10/10 02:33:34 jentron Exp $";
+static const char *IdSrc = "$Id: FGOutput.cpp,v 1.67 2012/04/27 12:14:56 jberndt Exp $";
static const char *IdHdr = ID_OUTPUT;
// (stolen from FGFS native_fdm.cxx)
}
if (SubSystems & ssMoments) {
outstream << delimeter;
- outstream << "L_{Aero} (ft-lbs)" + delimeter + "M_{Aero} ( ft-lbs)" + delimeter + "N_{Aero} (ft-lbs)" + delimeter;
+ outstream << "L_{Aero} (ft-lbs)" + delimeter + "M_{Aero} (ft-lbs)" + delimeter + "N_{Aero} (ft-lbs)" + delimeter;
outstream << "L_{Prop} (ft-lbs)" + delimeter + "M_{Prop} (ft-lbs)" + delimeter + "N_{Prop} (ft-lbs)" + delimeter;
outstream << "L_{Gear} (ft-lbs)" + delimeter + "M_{Gear} (ft-lbs)" + delimeter + "N_{Gear} (ft-lbs)" + delimeter;
outstream << "L_{ext} (ft-lbs)" + delimeter + "M_{ext} (ft-lbs)" + delimeter + "N_{ext} (ft-lbs)" + delimeter;
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+bool FGOutput::Load(int subSystems, std::string protocol, std::string type, std::string port, std::string name, double outRate, std::vector<FGPropertyManager *> & outputProperties)
+{
+ SetType(type);
+ SetRate(outRate);
+ SubSystems = subSystems;
+ OutputProperties = outputProperties;
+
+ if (((Type == otCSV) || (Type == otTab)) && (name != "cout") && (name !="COUT"))
+ name = FDMExec->GetRootDir() + name;
+
+ if (!port.empty() && (Type == otSocket || Type == otFlightGear)) {
+ SetProtocol(protocol);
+ socket = new FGfdmSocket(name, atoi(port.c_str()), Protocol);
+ } else {
+ BaseFilename = Filename = name;
+ }
+
+ Debug(2);
+ return true;
+}
+
bool FGOutput::Load(Element* element)
{
- string parameter="";
- string name="";
- double OutRate = 0.0;
- unsigned int port;
+ int subSystems = 0;
Element *property_element;
-
- string separator = "/";
+ std::vector<FGPropertyManager *> outputProperties;
if (!DirectivesFile.empty()) { // A directives filename from the command line overrides
output_file_name = DirectivesFile; // one found in the config file.
if (!document) return false;
- SetType(document->GetAttributeValue("type"));
-
- name = document->GetAttributeValue("name");
- if (((Type == otCSV) || (Type == otTab)) && (name != "cout") && (name !="COUT"))
- name = FDMExec->GetRootDir() + name;
-
- Port = document->GetAttributeValue("port");
- if (!Port.empty() && (Type == otSocket || Type == otFlightGear)) {
- port = atoi(Port.c_str());
- SetProtocol(document->GetAttributeValue("protocol"));
- socket = new FGfdmSocket(name, port, Protocol);
- } else {
- BaseFilename = Filename = name;
- }
+ string type = document->GetAttributeValue("type");
+ string name = document->GetAttributeValue("name");
+ string port = document->GetAttributeValue("port");
+ string protocol = document->GetAttributeValue("protocol");
if (!document->GetAttributeValue("rate").empty()) {
- OutRate = document->GetAttributeValueAsNumber("rate");
+ rate = document->GetAttributeValueAsNumber("rate");
} else {
- OutRate = 1;
+ rate = 1;
}
if (document->FindElementValue("simulation") == string("ON"))
- SubSystems += ssSimulation;
+ subSystems += ssSimulation;
if (document->FindElementValue("aerosurfaces") == string("ON"))
- SubSystems += ssAerosurfaces;
+ subSystems += ssAerosurfaces;
if (document->FindElementValue("rates") == string("ON"))
- SubSystems += ssRates;
+ subSystems += ssRates;
if (document->FindElementValue("velocities") == string("ON"))
- SubSystems += ssVelocities;
+ subSystems += ssVelocities;
if (document->FindElementValue("forces") == string("ON"))
- SubSystems += ssForces;
+ subSystems += ssForces;
if (document->FindElementValue("moments") == string("ON"))
- SubSystems += ssMoments;
+ subSystems += ssMoments;
if (document->FindElementValue("atmosphere") == string("ON"))
- SubSystems += ssAtmosphere;
+ subSystems += ssAtmosphere;
if (document->FindElementValue("massprops") == string("ON"))
- SubSystems += ssMassProps;
+ subSystems += ssMassProps;
if (document->FindElementValue("position") == string("ON"))
- SubSystems += ssPropagate;
- if (document->FindElementValue("coefficients") == string("ON"))
- SubSystems += ssAeroFunctions;
+ subSystems += ssPropagate;
+ if (document->FindElementValue("coefficients") == string("ON") || document->FindElementValue("aerodynamics") == string("ON"))
+ subSystems += ssAeroFunctions;
if (document->FindElementValue("ground_reactions") == string("ON"))
- SubSystems += ssGroundReactions;
+ subSystems += ssGroundReactions;
if (document->FindElementValue("fcs") == string("ON"))
- SubSystems += ssFCS;
+ subSystems += ssFCS;
if (document->FindElementValue("propulsion") == string("ON"))
- SubSystems += ssPropulsion;
+ subSystems += ssPropulsion;
property_element = document->FindElement("property");
while (property_element) {
string property_str = property_element->GetDataLine();
<< " not be logged. You should check your configuration file."
<< reset << endl;
} else {
- OutputProperties.push_back(node);
+ outputProperties.push_back(node);
}
property_element = document->FindNextElement("property");
}
- SetRate(OutRate);
+ return Load(subSystems, protocol, type, port, name, rate, outputProperties);
+}
- Debug(2);
- return true;
-}
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
DEFINITIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#define ID_OUTPUT "$Id: FGOutput.h,v 1.24 2011/11/10 12:06:14 jberndt Exp $"
+#define ID_OUTPUT "$Id: FGOutput.h,v 1.25 2012/02/07 23:15:37 bcoconni Exp $"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FORWARD DECLARATIONS
propulsion ON|OFF
</pre>
NOTE that Time is always output with the data.
- @version $Id: FGOutput.h,v 1.24 2011/11/10 12:06:14 jberndt Exp $
+ @version $Id: FGOutput.h,v 1.25 2012/02/07 23:15:37 bcoconni Exp $
*/
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
bool Toggle(void) {enabled = !enabled; return enabled;}
bool Load(Element* el);
+ bool Load(int subSystems, std::string protocol, std::string type, std::string port,
+ std::string name, double outRate,
+ std::vector<FGPropertyManager *> & outputProperties);
string GetOutputFileName(void) const {return Filename;}
/// Subsystem types for specifying which will be output in the FDM data logging
Wiley & Sons, 1979 ISBN 0-471-03032-5
[5] Bernard Etkin, "Dynamics of Flight, Stability and Control", Wiley & Sons,
1982 ISBN 0-471-08936-2
-[6] Erin Catto, "Iterative Dynamics with Temporal Coherence", February 22, 2005
+[6] S. Buss, "Accurate and Efficient Simulation of Rigid Body Rotations",
+ Technical Report, Department of Mathematics, University of California,
+ San Diego, 1999
+[7] Barker L.E., Bowles R.L. and Williams L.H., "Development and Application of
+ a Local Linearization Algorithm for the Integration of Quaternion Rate
+ Equations in Real-Time Flight Simulation Problems", NASA TN D-7347,
+ December 1973
+[8] Phillips W.F, Hailey C.E and Gebert G.A, "Review of Attitude Representations
+ Used for Aircraft Kinematics", Journal Of Aircraft Vol. 38, No. 4,
+ July-August 2001
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
INCLUDES
namespace JSBSim {
-static const char *IdSrc = "$Id: FGPropagate.cpp,v 1.100 2011/11/06 18:14:51 bcoconni Exp $";
+static const char *IdSrc = "$Id: FGPropagate.cpp,v 1.105 2012/03/26 21:26:11 bcoconni Exp $";
static const char *IdHdr = ID_PROPAGATE;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
{
Debug(0);
Name = "FGPropagate";
-
+
vInertialVelocity.InitMatrix();
/// These define the indices use to select the various integrators.
and current derivatives. Based on these values we compute an approximation to the
state values for (now + dt).
-In the code below, variables named beginning with a small "v" refer to a
+In the code below, variables named beginning with a small "v" refer to a
a column vector, variables named beginning with a "T" refer to a transformation
matrix. ECEF refers to Earth Centered Earth Fixed. ECI refers to Earth Centered
Inertial.
double dt = in.DeltaT * rate; // The 'stepsize'
- RunPreFunctions();
-
// Propagate rotational / translational velocity, angular /translational position, respectively.
- Integrate(VState.vPQRi, in.vPQRidot, VState.dqPQRidot, dt, integrator_rotational_rate);
Integrate(VState.qAttitudeECI, in.vQtrndot, VState.dqQtrndot, dt, integrator_rotational_position);
+ Integrate(VState.vPQRi, in.vPQRidot, VState.dqPQRidot, dt, integrator_rotational_rate);
Integrate(VState.vInertialPosition, VState.vInertialVelocity, VState.dqInertialVelocity, dt, integrator_translational_position);
Integrate(VState.vInertialVelocity, in.vUVWidot, VState.dqUVWidot, dt, integrator_translational_rate);
// vLocation vector.
UpdateLocationMatrices();
- // 5. Normalize the ECI Attitude quaternion
- VState.qAttitudeECI.Normalize();
-
- // 6. Update the "Orientation-based" transformation matrices from the updated
+ // 5. Update the "Orientation-based" transformation matrices from the updated
// orientation quaternion and vLocation vector.
UpdateBodyMatrices();
// Compute vehicle velocity wrt ECEF frame, expressed in Local horizontal frame.
vVel = Tb2l * VState.vUVW;
- RunPostFunctions();
-
Debug(2);
return false;
}
// Transform the velocity vector of the body relative to the origin (Earth
// center) to be expressed in the inertial frame, and add the vehicle velocity
// contribution due to the rotation of the planet.
- // Reference: See Stevens and Lewis, "Aircraft Control and Simulation",
+ // Reference: See Stevens and Lewis, "Aircraft Control and Simulation",
// Second edition (2004), eqn 1.5-16c (page 50)
void FGPropagate::CalculateInertialVelocity(void)
break;
case eNone: // do nothing, freeze translational rate
break;
+ default:
+ break;
}
}
break;
case eAdamsBashforth4: Integrand += (1/24.0)*dt*(55.0*ValDot[0] - 59.0*ValDot[1] + 37.0*ValDot[2] - 9.0*ValDot[3]);
break;
+ case eBuss1:
+ {
+ // This is the first order method as described in Samuel R. Buss paper[6].
+ // The formula from Buss' paper is transposed below to quaternions and is
+ // actually the exact solution of the quaternion differential equation
+ // qdot = 1/2*w*q when w is constant.
+ Integrand = Integrand * QExp(0.5 * dt * VState.vPQRi);
+ }
+ return; // No need to normalize since the quaternion exponential is always normal
+ case eBuss2:
+ {
+ // This is the 'augmented second-order method' from S.R. Buss paper [6].
+ // Unlike Runge-Kutta or Adams-Bashforth, it is a one-pass second-order
+ // method (see reference [6]).
+ FGColumnVector3 wi = VState.vPQRi;
+ FGColumnVector3 wdoti = in.vPQRidot;
+ FGColumnVector3 omega = wi + 0.5*dt*wdoti + dt*dt/12.*wdoti*wi;
+ Integrand = Integrand * QExp(0.5 * dt * omega);
+ }
+ return; // No need to normalize since the quaternion exponential is always normal
+ case eLocalLinearization:
+ {
+ // This is the local linearization algorithm of Barker et al. (see ref. [7])
+ // It is also a one-pass second-order method. The code below is based on the
+ // more compact formulation issued from equation (107) of ref. [8]. The
+ // constants C1, C2, C3 and C4 have the same value than those in ref. [7] pp. 11
+ FGColumnVector3 wi = 0.5 * VState.vPQRi;
+ FGColumnVector3 wdoti = 0.5 * in.vPQRidot;
+ double omegak2 = DotProduct(VState.vPQRi, VState.vPQRi);
+ double omegak = omegak2 > 1E-6 ? sqrt(omegak2) : 1E-6;
+ double rhok = 0.5 * dt * omegak;
+ double C1 = cos(rhok);
+ double C2 = 2.0 * sin(rhok) / omegak;
+ double C3 = 4.0 * (1.0 - C1) / (omegak*omegak);
+ double C4 = 4.0 * (dt - C2) / (omegak*omegak);
+ FGColumnVector3 Omega = C2*wi + C3*wdoti + C4*wi*wdoti;
+ FGQuaternion q;
+
+ q(1) = C1 - C4*DotProduct(wi, wdoti);
+ q(2) = Omega(eP);
+ q(3) = Omega(eQ);
+ q(4) = Omega(eR);
+
+ Integrand = Integrand * q;
+
+ /* Cross check with ref. [7] pp.11-12 formulas and code pp. 20
+ double pk = VState.vPQRi(eP);
+ double qk = VState.vPQRi(eQ);
+ double rk = VState.vPQRi(eR);
+ double pdotk = in.vPQRidot(eP);
+ double qdotk = in.vPQRidot(eQ);
+ double rdotk = in.vPQRidot(eR);
+ double Ap = -0.25 * (pk*pdotk + qk*qdotk + rk*rdotk);
+ double Bp = 0.25 * (pk*qdotk - qk*pdotk);
+ double Cp = 0.25 * (pdotk*rk - pk*rdotk);
+ double Dp = 0.25 * (qk*rdotk - qdotk*rk);
+ double C2p = sin(rhok) / omegak;
+ double C3p = 2.0 * (1.0 - cos(rhok)) / (omegak*omegak);
+ double H = C1 + C4 * Ap;
+ double G = -C2p*rk - C3p*rdotk + C4*Bp;
+ double J = C2p*qk + C3p*qdotk - C4*Cp;
+ double K = C2p*pk + C3p*pdotk - C4*Dp;
+
+ cout << "q: " << q << endl;
+
+ // Warning! In the paper of Barker et al. the quaternion components are not
+ // ordered the same way as in JSBSim (see equations (2) and (3) of ref. [7]
+ // as well as the comment just below equation (3))
+ cout << "FORTRAN: " << H << " , " << K << " , " << J << " , " << -G << endl;*/
+ }
+ break; // The quaternion q is not normal so the normalization needs to be done.
case eNone: // do nothing, freeze rotational rate
break;
+ default:
+ break;
}
+
+ Integrand.Normalize();
}
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
{
//ToDo: Shouldn't all of these be set from the vstate vector passed in?
VState.vLocation = vstate.vLocation;
- VState.vLocation.SetEarthPositionAngle(vstate.vLocation.GetEPA());
Ti2ec = VState.vLocation.GetTi2ec(); // useless ?
Tec2i = Ti2ec.Transposed();
UpdateLocationMatrices();
PropertyManager->Tie("attitude/roll-rad", this, (int)ePhi, (PMF)&FGPropagate::GetEuler);
PropertyManager->Tie("attitude/pitch-rad", this, (int)eTht, (PMF)&FGPropagate::GetEuler);
PropertyManager->Tie("attitude/heading-true-rad", this, (int)ePsi, (PMF)&FGPropagate::GetEuler);
-
+
PropertyManager->Tie("simulation/integrator/rate/rotational", (int*)&integrator_rotational_rate);
PropertyManager->Tie("simulation/integrator/rate/translational", (int*)&integrator_translational_rate);
PropertyManager->Tie("simulation/integrator/position/rotational", (int*)&integrator_rotational_position);
if (debug_lvl & 4 ) { // Run() method entry print for FGModel-derived objects
}
if (debug_lvl & 8 && from == 2) { // Runtime state variables
- cout << endl << fgblue << highint << left
+ cout << endl << fgblue << highint << left
<< " Propagation Report (English units: ft, degrees) at simulation time " << FDMExec->GetSimTime() << " seconds"
<< reset << endl;
cout << endl;
DEFINITIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#define ID_PROPAGATE "$Id: FGPropagate.h,v 1.67 2011/11/09 22:07:17 bcoconni Exp $"
+#define ID_PROPAGATE "$Id: FGPropagate.h,v 1.69 2012/04/29 13:27:51 bcoconni Exp $"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FORWARD DECLARATIONS
state of the vehicle given the forces and moments that act on it. The
integration accounts for a rotating Earth.
- Integration of rotational and translation position and rate can be
+ Integration of rotational and translation position and rate can be
customized as needed or frozen by the selection of no integrator. The
- selection of which integrator to use is done through the setting of
+ selection of which integrator to use is done through the setting of
the associated property. There are four properties which can be set:
-
+
@code
simulation/integrator/rate/rotational
simulation/integrator/rate/translational
simulation/integrator/position/rotational
simulation/integrator/position/translational
@endcode
-
+
Each of the integrators listed above can be set to one of the following values:
@code
@endcode
@author Jon S. Berndt, Mathias Froehlich, Bertrand Coconnier
- @version $Id: FGPropagate.h,v 1.67 2011/11/09 22:07:17 bcoconni Exp $
+ @version $Id: FGPropagate.h,v 1.69 2012/04/29 13:27:51 bcoconni Exp $
*/
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
/// Destructor
~FGPropagate();
-
+
/// These define the indices use to select the various integrators.
- enum eIntegrateType {eNone = 0, eRectEuler, eTrapezoidal, eAdamsBashforth2, eAdamsBashforth3, eAdamsBashforth4};
+ enum eIntegrateType {eNone = 0, eRectEuler, eTrapezoidal, eAdamsBashforth2,
+ eAdamsBashforth3, eAdamsBashforth4, eBuss1, eBuss2, eLocalLinearization};
/** Initializes the FGPropagate class after instantiation and prior to first execution.
- The base class FGModel::InitModel is called first, initializing pointers to the
+ The base class FGModel::InitModel is called first, initializing pointers to the
other FGModel objects (and others). */
bool InitModel(void);
/** Runs the state propagation model; called by the Executive
Can pass in a value indicating if the executive is directing the simulation to Hold.
- @param Holding if true, the executive has been directed to hold the sim from
+ @param Holding if true, the executive has been directed to hold the sim from
advancing time. Some models may ignore this flag, such as the Input
model, which may need to be active to listen on a socket for the
"Resume" command to be given.
expressed in Local horizontal frame.
*/
const FGColumnVector3& GetVel(void) const { return vVel; }
-
+
/** Retrieves the body frame vehicle velocity vector.
The vector returned is represented by an FGColumnVector reference. The vector
for the velocity in Body frame is organized (Vx, Vy, Vz). The vector
@return The body frame vehicle velocity vector in ft/sec.
*/
const FGColumnVector3& GetUVW(void) const { return VState.vUVW; }
-
+
/** Retrieves the body angular rates vector, relative to the ECEF frame.
Retrieves the body angular rates (p, q, r), which are calculated by integration
of the angular acceleration.
@return The body frame angular rates in rad/sec.
*/
const FGColumnVector3& GetPQR(void) const {return VState.vPQR;}
-
+
/** Retrieves the body angular rates vector, relative to the ECI (inertial) frame.
Retrieves the body angular rates (p, q, r), which are calculated by integration
of the angular acceleration.
FGColumnVector3 vVel;
FGColumnVector3 vInertialVelocity;
FGColumnVector3 vLocation;
- FGColumnVector3 vDeltaXYZEC;
FGMatrix33 Tec2b;
FGMatrix33 Tb2ec;
FGMatrix33 Tl2b; // local to body frame matrix copy for immediate local use
FGMatrix33 Tb2i; // body to ECI frame rotation matrix
FGMatrix33 Ti2l;
FGMatrix33 Tl2i;
-
+
double VehicleRadius;
FGColumnVector3 LocalTerrainVelocity, LocalTerrainAngularVelocity;
namespace JSBSim {
-static const char *IdSrc = "$Id: FGPropulsion.cpp,v 1.52 2011/10/31 14:54:41 bcoconni Exp $";
+static const char *IdSrc = "$Id: FGPropulsion.cpp,v 1.61 2012/04/14 18:10:44 bcoconni Exp $";
static const char *IdHdr = ID_PROPULSION;
extern short debug_lvl;
unsigned int TanksWithOxidizer=0, CurrentOxidizerTankPriority=1;
vector <int> FeedListFuel, FeedListOxi;
bool Starved = true; // Initially set Starved to true. Set to false in code below.
-// bool hasOxTanks = false;
+ bool hasOxTanks = false;
// For this engine,
// 1) Count how many fuel tanks with the current priority level have fuel
if (TanksWithFuel == 0) CurrentFuelTankPriority++; // No tanks at this priority, try next priority
}
+ bool FuelStarved = Starved;
+ Starved = true;
+
// Process Oxidizer tanks, if any
if (engine->GetType() == FGEngine::etRocket) {
while ((TanksWithOxidizer == 0) && (CurrentOxidizerTankPriority <= numTanks)) {
// Skip this here (done above)
break;
case FGTank::ttOXIDIZER:
-// hasOxTanks = true;
+ hasOxTanks = true;
if (Tank->GetContents() > 0.0 && Tank->GetSelected() && TankPriority == CurrentOxidizerTankPriority) {
TanksWithOxidizer++;
if (TanksWithFuel > 0) Starved = false;
}
}
- engine->SetStarved(Starved); // Tanks can be refilled, so be sure to reset engine Starved flag here.
+ bool OxiStarved = Starved;
+
+ engine->SetStarved(FuelStarved || (hasOxTanks && OxiStarved)); // Tanks can be refilled, so be sure to reset engine Starved flag here.
// No fuel or fuel/oxidizer found at any priority!
- if (Starved) return;
+// if (Starved) return;
+ if (FuelStarved || (hasOxTanks && OxiStarved)) return;
double FuelToBurn = engine->CalcFuelNeed(); // How much fuel does this engine need?
double FuelNeededPerTank = FuelToBurn / TanksWithFuel; // Determine fuel needed per tank.
fullpath = enginePath + separator;
localpath = aircraftPath + separator + "Engines" + separator;
- engine_file.open(string(fullpath + engine_filename + ".xml").c_str());
+ engine_file.open(string(localpath + engine_filename + ".xml").c_str());
if ( !engine_file.is_open()) {
- engine_file.open(string(localpath + engine_filename + ".xml").c_str());
+ engine_file.open(string(fullpath + engine_filename + ".xml").c_str());
if ( !engine_file.is_open()) {
cerr << " Could not open engine file: " << engine_filename << " in path "
<< fullpath << " or " << localpath << endl;
return string("");
} else {
- return string(localpath + engine_filename + ".xml");
+ return string(fullpath + engine_filename + ".xml");
}
}
- return string(fullpath + engine_filename + ".xml");
+ return string(localpath + engine_filename + ".xml");
}
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
fullpath = enginePath + separator;
localpath = aircraftPath + separator + "Engines" + separator;
- engine_file->open(string(fullpath + engine_filename + ".xml").c_str());
+ engine_file->open(string(localpath + engine_filename + ".xml").c_str());
if ( !engine_file->is_open()) {
- engine_file->open(string(localpath + engine_filename + ".xml").c_str());
+ engine_file->open(string(fullpath + engine_filename + ".xml").c_str());
if ( !engine_file->is_open()) {
cerr << " Could not open engine file: " << engine_filename << " in path "
- << fullpath << " or " << localpath << endl;
+ << localpath << " or " << fullpath << endl;
}
}
return engine_file;
else if (Tanks[i]->GetType() == FGTank::ttOXIDIZER) buf << delimiter << "Oxidizer Tank " << i;
}
+ PropulsionStrings += buf.str();
+ buf.str("");
+
return PropulsionStrings;
}
buf << Tanks[i]->GetContents();
}
+ PropulsionValues += buf.str();
+ buf.str("");
+
return PropulsionValues;
}
{
vXYZtank_arm.InitMatrix();
for (unsigned int i=0; i<Tanks.size(); i++) {
- vXYZtank_arm(eX) += Tanks[i]->GetXYZ(eX) * Tanks[i]->GetContents();
- vXYZtank_arm(eY) += Tanks[i]->GetXYZ(eY) * Tanks[i]->GetContents();
- vXYZtank_arm(eZ) += Tanks[i]->GetXYZ(eZ) * Tanks[i]->GetContents();
+ vXYZtank_arm += Tanks[i]->GetXYZ() * Tanks[i]->GetContents();
}
return vXYZtank_arm;
}
const FGMatrix33& FGPropulsion::CalculateTankInertias(void)
{
+ const FGMatrix33 Ts2b(-inchtoft, 0., 0., 0., inchtoft, 0., 0., 0., -inchtoft);
unsigned int size;
size = Tanks.size();
tankJ = FGMatrix33();
for (unsigned int i=0; i<size; i++) {
+ FGColumnVector3 vTankBodyVec = Ts2b * (in.vXYZcg - Tanks[i]->GetXYZ());
+
tankJ += FDMExec->GetMassBalance()->GetPointmassInertia( lbtoslug * Tanks[i]->GetContents(),
- Tanks[i]->GetXYZ() );
+ vTankBodyVec);
tankJ(1,1) += Tanks[i]->GetIxx();
tankJ(2,2) += Tanks[i]->GetIyy();
tankJ(3,3) += Tanks[i]->GetIzz();
Module: FGMSIS.cpp
Author: David Culp
- (incorporated into C++ JSBSim class heirarchy, see model authors below)
+ (incorporated into C++ JSBSim class hierarchy, see model authors below)
Date started: 12/14/03
Purpose: Models the MSIS-00 atmosphere
namespace JSBSim {
-static const char *IdSrc = "$Id: FGMSIS.cpp,v 1.18 2011/06/21 13:54:40 jberndt Exp $";
+static const char *IdSrc = "$Id: FGMSIS.cpp,v 1.19 2011/12/11 17:03:05 bcoconni Exp $";
static const char *IdHdr = ID_MSIS;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if (FGModel::Run(Holding)) return true;
if (Holding) return false;
- RunPreFunctions();
-
double h = FDMExec->GetPropagate()->GetAltitudeASL();
//do temp, pressure, and density first
-// if (!useExternal) {
+ //if (!useExternal) {
// get sea-level values
Calculate(FDMExec->GetAuxiliary()->GetDayOfYear(),
FDMExec->GetAuxiliary()->GetSecondsInDay(),
h,
FDMExec->GetPropagate()->GetLocation().GetLatitudeDeg(),
FDMExec->GetPropagate()->GetLocation().GetLongitudeDeg());
-// intTemperature = output.t[1] * 1.8;
-// intDensity = output.d[5] * 1.940321;
-// intPressure = 1716.488 * intDensity * intTemperature;
+ //intTemperature = output.t[1] * 1.8;
+ //intDensity = output.d[5] * 1.940321;
+ //intPressure = 1716.488 * intDensity * intTemperature;
//cout << "T=" << intTemperature << " D=" << intDensity << " P=";
//cout << intPressure << " a=" << soundspeed << endl;
-// }
-
- RunPostFunctions();
+ //}
Debug(2);
namespace JSBSim {
-static const char *IdSrc = "$Id: FGStandardAtmosphere.cpp,v 1.20 2011/09/18 12:06:21 bcoconni Exp $";
+static const char *IdSrc = "$Id: FGStandardAtmosphere.cpp,v 1.21 2012/04/13 13:18:27 jberndt Exp $";
static const char *IdHdr = ID_STANDARDATMOSPHERE;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-void FGStandardAtmosphere::SetPressureSL(double pressure, ePressure unit)
+void FGStandardAtmosphere::SetPressureSL(ePressure unit, double pressure)
{
double press = ConvertToPSF(pressure, unit);
PropertyManager->Tie("atmosphere/SL-graded-delta-T", this, eRankine,
(PMFi)&FGStandardAtmosphere::GetTemperatureDeltaGradient,
(PMF)&FGStandardAtmosphere::SetSLTemperatureGradedDelta);
+ PropertyManager->Tie("atmosphere/P-sl-psf", this, ePSF,
+ (PMFi)&FGStandardAtmosphere::GetPressureSL,
+ (PMF)&FGStandardAtmosphere::SetPressureSL);
}
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
DEFINITIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#define ID_STANDARDATMOSPHERE "$Id: FGStandardAtmosphere.h,v 1.16 2011/09/18 12:06:21 bcoconni Exp $"
+#define ID_STANDARDATMOSPHERE "$Id: FGStandardAtmosphere.h,v 1.17 2012/04/13 13:18:27 jberndt Exp $"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FORWARD DECLARATIONS
@author Jon Berndt
@see "U.S. Standard Atmosphere, 1976", NASA TM-X-74335
- @version $Id: FGStandardAtmosphere.h,v 1.16 2011/09/18 12:06:21 bcoconni Exp $
+ @version $Id: FGStandardAtmosphere.h,v 1.17 2012/04/13 13:18:27 jberndt Exp $
*/
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
/** Sets the sea level pressure for modeling an off-standard pressure
profile. This could be useful in the case where the pressure at an
airfield is known or set for a particular simulation run.
- @param pressure The pressure in the units specified (PSF by default).
+ @param pressure The pressure in the units specified.
@param unit the unit of measure that the specified pressure is
supplied in.*/
- virtual void SetPressureSL(double pressure, ePressure unit=ePSF);
+ virtual void SetPressureSL(ePressure unit, double pressure);
/** Resets the sea level to the Standard sea level pressure, and recalculates
dependent parameters so that the pressure calculations are standard. */
namespace JSBSim {
-static const char *IdSrc = "$Id: FGWinds.cpp,v 1.6 2011/11/10 12:02:34 jberndt Exp $";
+static const char *IdSrc = "$Id: FGWinds.cpp,v 1.7 2011/12/11 17:03:05 bcoconni Exp $";
static const char *IdHdr = ID_WINDS;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if (FGModel::Run(Holding)) return true;
if (Holding) return false;
- RunPreFunctions();
-
if (turbType != ttNone) Turbulence(in.AltitudeASL);
if (oneMinusCosineGust.gustProfile.Running) CosineGust();
if (vWindNED(eX) != 0.0) psiw = atan2( vWindNED(eY), vWindNED(eX) );
if (psiw < 0) psiw += 2*M_PI;
- RunPostFunctions();
-
Debug(2);
return false;
}
{
double mag = GetWindspeed();
psiw = dir;
- SetWindspeed(mag);
+ SetWindspeed(mag);
}
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
{
switch (turbType) {
- case ttCulp: {
+ case ttCulp: {
vTurbPQR(eP) = wind_from_clockwise;
if (TurbGain == 0.0) return;
-
+
// keep the inputs within allowable limts for this model
if (TurbGain < 0.0) TurbGain = 0.0;
if (TurbGain > 1.0) TurbGain = 1.0;
if (time > target_time) {
spike = 1.0;
target_time = 0.0;
- }
+ }
// max vertical wind speed in fps, corresponds to TurbGain = 1.0
double max_vs = 40;
vTurbulenceNED(3)+= delta;
if (in.DistanceAGL/in.wingspan < 3.0)
vTurbulenceNED(3) *= in.DistanceAGL/in.wingspan * 0.3333;
-
+
// Yaw component of turbulence.
vTurbulenceNED(1) = sin( delta * 3.0 );
vTurbulenceNED(2) = cos( delta * 3.0 );
(PMFd)&FGWinds::SetGustNED);
PropertyManager->Tie("atmosphere/gust-down-fps", this, eDown, (PMF)&FGWinds::GetGustNED,
(PMFd)&FGWinds::SetGustNED);
-
+
// User-specified 1 - cosine gust parameters (in specified frame)
PropertyManager->Tie("atmosphere/cosine-gust/startup-duration-sec", this, (Ptr)0L, &FGWinds::StartupGustDuration);
PropertyManager->Tie("atmosphere/cosine-gust/steady-duration-sec", this, (Ptr)0L, &FGWinds::SteadyGustDuration);
}
if (debug_lvl & 16) { // Sanity checking
}
- if (debug_lvl & 128) { //
+ if (debug_lvl & 128) { //
}
if (debug_lvl & 64) {
if (from == 0) { // Constructor
DEFINITIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#define ID_ACCELEROMETER "$Id: FGAccelerometer.h,v 1.5 2011/07/17 13:51:23 jberndt Exp $"
+#define ID_ACCELEROMETER "$Id: FGAccelerometer.h,v 1.6 2012/01/08 12:39:14 bcoconni Exp $"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FORWARD DECLARATIONS
@code
<accelerometer name="name">
- <input> property </input>
+ <location unit="{IN | M}">
+ <x> number </x>
+ <y> number </y>
+ <z> number </z>
+ </location>
+ <orientation unit="{RAD | DEG}">
+ <pitch> {number} </pitch>
+ <roll> {number} </roll>
+ <yaw> {number} </yaw>
+ </orientation>
+ <axis> {X | Y | Z} </axis>
<lag> number </lag>
<noise variation="PERCENT|ABSOLUTE"> number </noise>
<quantization name="name">
<max> number </max>
</quantization>
<drift_rate> number </drift_rate>
+ <gain> number </gain>
<bias> number </bias>
</accelerometer>
@endcode
Example:
@code
-<accelerometer name="aero/accelerometer/qbar">
- <input> aero/qbar </input>
+<accelerometer name="aero/accelerometer/right_tip_wing">
+ <location unit="IN">
+ <x> 43.2 </x>
+ <y> 214. </y>
+ <z> 59.4 </z>
+ </location>
+ <axis> Z </axis>
<lag> 0.5 </lag>
<noise variation="PERCENT"> 2 </noise>
- <quantization name="aero/accelerometer/quantized/qbar">
+ <quantization name="aero/accelerometer/quantized/right_tip_wing">
<bits> 12 </bits>
<min> 0 </min>
<max> 400 </max>
time.
@author Jon S. Berndt
-@version $Revision: 1.5 $
+@version $Revision: 1.6 $
*/
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
namespace JSBSim {
-static const char *IdSrc = "$Id: FGActuator.cpp,v 1.22 2011/07/12 21:40:32 jentron Exp $";
+static const char *IdSrc = "$Id: FGActuator.cpp,v 1.23 2012/04/08 15:04:41 jberndt Exp $";
static const char *IdHdr = ID_ACTUATOR;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
fail_zero = fail_hardover = fail_stuck = false;
ca = cb = 0.0;
initialized = 0;
+ saturated = false;
if ( element->FindElement("deadband_width") ) {
deadband_width = element->FindElementValueAsNumber("deadband_width");
initialized = 1;
Clip();
+
+ if (clip) {
+ saturated = false;
+ if (Output >= clipmax) saturated = true;
+ else if (Output <= clipmin) saturated = true;
+ }
+
if (IsOutput) SetOutput();
return true;
const string tmp_zero = tmp + "/malfunction/fail_zero";
const string tmp_hardover = tmp + "/malfunction/fail_hardover";
const string tmp_stuck = tmp + "/malfunction/fail_stuck";
+ const string tmp_sat = tmp + "/saturated";
PropertyManager->Tie( tmp_zero, this, &FGActuator::GetFailZero, &FGActuator::SetFailZero);
PropertyManager->Tie( tmp_hardover, this, &FGActuator::GetFailHardover, &FGActuator::SetFailHardover);
PropertyManager->Tie( tmp_stuck, this, &FGActuator::GetFailStuck, &FGActuator::SetFailStuck);
+ PropertyManager->Tie( tmp_sat, this, &FGActuator::IsSaturated);
}
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
DEFINITIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#define ID_ACTUATOR "$Id: FGActuator.h,v 1.12 2011/07/12 21:40:32 jentron Exp $"
+#define ID_ACTUATOR "$Id: FGActuator.h,v 1.13 2012/04/08 15:04:41 jberndt Exp $"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FORWARD DECLARATIONS
@endcode
@author Jon S. Berndt
-@version $Revision: 1.12 $
+@version $Revision: 1.13 $
*/
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
/** This function fails the actuator to zero. The motion to zero
will flow through the lag, hysteresis, and rate limiting
functions if those are activated. */
- inline void SetFailZero(bool set) {fail_zero = set;}
- inline void SetFailHardover(bool set) {fail_hardover = set;}
- inline void SetFailStuck(bool set) {fail_stuck = set;}
-
- inline bool GetFailZero(void) const {return fail_zero;}
- inline bool GetFailHardover(void) const {return fail_hardover;}
- inline bool GetFailStuck(void) const {return fail_stuck;}
+ void SetFailZero(bool set) {fail_zero = set;}
+ void SetFailHardover(bool set) {fail_hardover = set;}
+ void SetFailStuck(bool set) {fail_stuck = set;}
+
+ bool GetFailZero(void) const {return fail_zero;}
+ bool GetFailHardover(void) const {return fail_hardover;}
+ bool GetFailStuck(void) const {return fail_stuck;}
+ bool IsSaturated(void) const {return saturated;}
private:
double span;
bool fail_hardover;
bool fail_stuck;
bool initialized;
+ bool saturated;
void Hysteresis(void);
void Lag(void);
namespace JSBSim {
-static const char *IdSrc = "$Id: FGPID.cpp,v 1.19 2011/05/05 11:44:11 jberndt Exp $";
+static const char *IdSrc = "$Id: FGPID.cpp,v 1.20 2012/05/10 12:10:48 jberndt Exp $";
static const char *IdHdr = ID_PID;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
I_out_total = 0.0;
Input_prev = Input_prev2 = 0.0;
Trigger = 0;
+ ProcessVariableDot = 0;
+ IsStandard = false;
+ IntType = eNone; // No integrator initially defined.
+
+ string pid_type = element->GetAttributeValue("type");
+
+ if (pid_type == "standard") IsStandard = true;
if ( element->FindElement("kp") ) {
kp_string = element->FindElementValue("kp");
if ( element->FindElement("ki") ) {
ki_string = element->FindElementValue("ki");
+
+ string integ_type = element->FindElement("ki")->GetAttributeValue("type");
+ if (integ_type == "rect") { // Use rectangular integration
+ IntType = eRectEuler;
+ } else if (integ_type == "trap") { // Use trapezoidal integration
+ IntType = eTrapezoidal;
+ } else if (integ_type == "ab2") { // Use Adams Bashforth 2nd order integration
+ IntType = eAdamsBashforth2;
+ } else if (integ_type == "ab3") { // Use Adams Bashforth 3rd order integration
+ IntType = eAdamsBashforth3;
+ } else { // Use default Adams Bashforth 2nd order integration
+ IntType = eAdamsBashforth2;
+ }
+
if (!is_number(ki_string)) { // property
if (ki_string[0] == '-') {
KiPropertySign = -1.0;
}
}
+ if (element->FindElement("pvdot")) {
+ ProcessVariableDot = PropertyManager->GetNode(element->FindElementValue("pvdot"));
+ }
+
if (element->FindElement("trigger")) {
Trigger = PropertyManager->GetNode(element->FindElementValue("trigger"));
}
bool FGPID::Run(void )
{
double I_out_delta = 0.0;
- double P_out, D_out;
+ double Dval = 0;
Input = InputNodes[0]->getDoubleValue() * InputSigns[0];
if (KiPropertyNode != 0) Ki = KiPropertyNode->getDoubleValue() * KiPropertySign;
if (KdPropertyNode != 0) Kd = KdPropertyNode->getDoubleValue() * KdPropertySign;
- P_out = Kp * Input;
- D_out = (Kd / dt) * (Input - Input_prev);
+ if (ProcessVariableDot) {
+ Dval = ProcessVariableDot->getDoubleValue();
+ } else {
+ Dval = (Input - Input_prev)/dt;
+ }
// Do not continue to integrate the input to the integrator if a wind-up
// condition is sensed - that is, if the property pointed to by the trigger
// element is non-zero. Reset the integrator to 0.0 if the Trigger value
// is negative.
- if (Trigger != 0) {
- double test = Trigger->getDoubleValue();
- if (fabs(test) < 0.000001) {
- // I_out_delta = Ki * dt * Input; // Normal rectangular integrator
+ double test = 0.0;
+ if (Trigger != 0) test = Trigger->getDoubleValue();
+
+ if (fabs(test) < 0.000001) {
+ switch(IntType) {
+ case eRectEuler:
+ I_out_delta = Ki * dt * Input; // Normal rectangular integrator
+ break;
+ case eTrapezoidal:
+ I_out_delta = (Ki/2.0) * dt * (Input + Input_prev); // Trapezoidal integrator
+ break;
+ case eAdamsBashforth2:
I_out_delta = Ki * dt * (1.5*Input - 0.5*Input_prev); // 2nd order Adams Bashforth integrator
+ break;
+ case eAdamsBashforth3: // 3rd order Adams Bashforth integrator
+ I_out_delta = (Ki/12.0) * dt * (23.0*Input - 16.0*Input_prev + 5.0*Input_prev2);
+ break;
+ case eNone:
+ // No integator is defined or used.
+ I_out_delta = 0.0;
+ break;
}
- if (test < 0.0) I_out_total = 0.0; // Reset integrator to 0.0
- } else { // no anti-wind-up trigger defined
- // I_out_delta = Ki * dt * Input;
- I_out_delta = Ki * dt * (1.5*Input - 0.5*Input_prev); // 2nd order Adams Bashforth integrator
}
+
+ if (test < 0.0) I_out_total = 0.0; // Reset integrator to 0.0
I_out_total += I_out_delta;
- Output = P_out + I_out_total + D_out;
-
+ if (IsStandard) {
+ Output = Kp * (Input + I_out_total + Kd*Dval);
+ } else {
+ Output = Kp*Input + I_out_total + Kd*Dval;
+ }
+
Input_prev = Input;
Input_prev2 = Input_prev;
DEFINITIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#define ID_PID "$Id: FGPID.h,v 1.12 2009/10/24 22:59:30 jberndt Exp $"
+#define ID_PID "$Id: FGPID.h,v 1.13 2012/05/10 12:10:48 jberndt Exp $"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FORWARD DECLARATIONS
<h3>Configuration Format:</h3>
@code
-<pid name="{string}">
- <kp> {number} </kp>
- <ki> {number} </ki>
- <kd> {number} </kd>
- <trigger> {string} </trigger>
+<pid name="{string}" [type="standard"]>
+ <kp> {number|property} </kp>
+ <ki> {number|property} </ki>
+ <kd> {number|property} </kd>
+ <trigger> {property} </trigger>
</pid>
@endcode
+For the integration constant element, one can also supply the type attribute for
+what kind of integrator to be used, one of:
+
+- rect, for a rectangular integrator
+- trap, for a trapezoidal integrator
+- ab2, for a second order Adams Bashforth integrator
+- ab3, for a third order Adams Bashforth integrator
+
+For example,
+
+@code
+<pid name="fcs/heading-control">
+ <kp> 3 </kp>
+ <ki type="ab3"> 1 </ki>
+ <kd> 1 </kd>
+</pid>
+@endcode
+
+
+
<h3>Configuration Parameters:</h3>
<pre>
+ The values of kp, ki, and kd have slightly different interpretations depending on
+ whether the PID controller is a standard one, or an ideal/parallel one - with the latter
+ being the default.
+
kp - Proportional constant, default value 0.
ki - Integrative constant, default value 0.
kd - Derivative constant, default value 0.
- trigger - Property which is used to sense wind-up, optional.
+ trigger - Property which is used to sense wind-up, optional. Most often, the trigger
+ will be driven by the "saturated" property of a particular actuator. When
+ the relevant actuator has reached it's limits (if there are any, specified
+ by the <clipto> element) the automatically generated saturated property will
+ be greater than zero (true). If this property is used as the trigger for the
+ integrator, the integrator will not continue to integrate while the property
+ is still true (> 1), preventing wind-up.
+ pvdot - The property to be used as the process variable time derivative.
+
+
</pre>
@author Jon S. Berndt
- @version $Revision: 1.12 $
+ @version $Revision: 1.13 $
*/
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
bool Run (void);
void ResetPastStates(void) {Input_prev = Input_prev2 = Output = I_out_total = 0.0;}
+ /// These define the indices use to select the various integrators.
+ enum eIntegrateType {eNone = 0, eRectEuler, eTrapezoidal, eAdamsBashforth2, eAdamsBashforth3};
+
private:
- FGPropertyManager *Trigger;
double Kp, Ki, Kd;
double I_out_total;
double Input_prev, Input_prev2;
double KpPropertySign;
double KiPropertySign;
double KdPropertySign;
+
+ bool IsStandard;
+
+ eIntegrateType IntType;
+
+ FGPropertyManager *Trigger;
FGPropertyManager* KpPropertyNode;
FGPropertyManager* KiPropertyNode;
FGPropertyManager* KdPropertyNode;
+ FGPropertyManager* ProcessVariableDot;
void Debug(int from);
};
DEFINITIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#define ID_SENSOR "$Id: FGSensor.h,v 1.20 2011/08/18 12:42:17 jberndt Exp $"
+#define ID_SENSOR "$Id: FGSensor.h,v 1.21 2012/01/08 12:39:14 bcoconni Exp $"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FORWARD DECLARATIONS
<max> number </max>
</quantization>
<drift_rate> number </drift_rate>
+ <gain> number </gain>
<bias> number </bias>
<delay> number < /delay>
</sensor>
the number of frames to delay the output signal.
@author Jon S. Berndt
-@version $Revision: 1.20 $
+@version $Revision: 1.21 $
*/
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
namespace JSBSim {
-static const char *IdSrc = "$Id: FGEngine.cpp,v 1.48 2011/10/31 14:54:41 bcoconni Exp $";
+static const char *IdSrc = "$Id: FGEngine.cpp,v 1.50 2012/03/17 20:46:29 jentron Exp $";
static const char *IdHdr = ID_ENGINE;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
thruster_filename = thruster_element->GetAttributeValue("file");
if ( !thruster_filename.empty()) {
- thruster_fullpathname = fullpath + thruster_filename + ".xml";
+ thruster_fullpathname = localpath + thruster_filename + ".xml";
thruster_file.open(thruster_fullpathname.c_str());
if ( !thruster_file.is_open()) {
- thruster_fullpathname = localpath + thruster_filename + ".xml";
+ thruster_fullpathname = fullpath + thruster_filename + ".xml";
thruster_file.open(thruster_fullpathname.c_str());
if ( !thruster_file.is_open()) {
cerr << "Could not open thruster file: " << thruster_filename << ".xml" << endl;
DEFINITIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#define ID_ENGINE "$Id: FGEngine.h,v 1.29 2011/11/10 12:06:14 jberndt Exp $"
+#define ID_ENGINE "$Id: FGEngine.h,v 1.35 2012/04/14 18:10:44 bcoconni Exp $"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FORWARD DECLARATIONS
documentation for engine and thruster classes.
</pre>
@author Jon S. Berndt
- @version $Id: FGEngine.h,v 1.29 2011/11/10 12:06:14 jberndt Exp $
+ @version $Id: FGEngine.h,v 1.35 2012/04/14 18:10:44 bcoconni Exp $
*/
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
vector <double> MixturePos;
vector <double> PropAdvance;
vector <bool> PropFeather;
+ FGColumnVector3 vXYZcg; // CG coordinates expressed in the structural frame
double TotalDeltaT;
};
DEFINITIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#define ID_FORCE "$Id: FGForce.h,v 1.14 2011/10/31 14:54:41 bcoconni Exp $"
+#define ID_FORCE "$Id: FGForce.h,v 1.17 2012/04/01 17:05:51 bcoconni Exp $"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FORWARD DECLARATIONS
<br><br></p>
@author Tony Peden
- @version $Id: FGForce.h,v 1.14 2011/10/31 14:54:41 bcoconni Exp $
+ @version $Id: FGForce.h,v 1.17 2012/04/01 17:05:51 bcoconni Exp $
*/
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
TransformType ttype;
FGColumnVector3 vXYZn;
FGColumnVector3 vActingXYZn;
+ FGMatrix33 mT;
private:
FGColumnVector3 vFb;
FGColumnVector3 vM;
FGColumnVector3 vDXYZ;
- FGMatrix33 mT;
-
void Debug(int from);
};
}
namespace JSBSim {
-static const char *IdSrc = "$Id: FGNozzle.cpp,v 1.14 2011/08/03 03:21:06 jberndt Exp $";
+static const char *IdSrc = "$Id: FGNozzle.cpp,v 1.15 2012/03/18 15:48:35 jentron Exp $";
static const char *IdHdr = ID_NOZZLE;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-string FGNozzle::GetThrusterLabels(int id, string delimeter)
+string FGNozzle::GetThrusterLabels(int id, const string& delimeter)
{
std::ostringstream buf;
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-string FGNozzle::GetThrusterValues(int id, string delimeter)
+string FGNozzle::GetThrusterValues(int id, const string& delimeter)
{
std::ostringstream buf;
DEFINITIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#define ID_NOZZLE "$Id: FGNozzle.h,v 1.9 2011/08/03 03:21:06 jberndt Exp $";
+#define ID_NOZZLE "$Id: FGNozzle.h,v 1.10 2012/03/18 15:48:36 jentron Exp $";
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FORWARD DECLARATIONS
All parameters MUST be specified.
@author Jon S. Berndt
- @version $Id: FGNozzle.h,v 1.9 2011/08/03 03:21:06 jberndt Exp $
+ @version $Id: FGNozzle.h,v 1.10 2012/03/18 15:48:36 jentron Exp $
*/
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
~FGNozzle();
double Calculate(double vacThrust);
- string GetThrusterLabels(int id, string delimeter);
- string GetThrusterValues(int id, string delimeter);
+ string GetThrusterLabels(int id, const string& delimeter);
+ string GetThrusterValues(int id, const string& delimeter);
private:
// double PE;
namespace JSBSim {
-static const char *IdSrc = "$Id: FGPiston.cpp,v 1.68 2011/10/11 15:13:34 jentron Exp $";
+static const char *IdSrc = "$Id: FGPiston.cpp,v 1.71 2012/04/07 01:50:54 jentron Exp $";
static const char *IdHdr = ID_PISTON;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
MaxHP = 200;
MinManifoldPressure_inHg = 6.5;
MaxManifoldPressure_inHg = 28.5;
+ ManifoldPressureLag=1.0;
ISFC = -1;
volumetric_efficiency = 0.85;
Bore = 5.125;
FMEPStatic = 46500;
Cooling_Factor = 0.5144444;
StaticFriction_HP = 1.5;
+ StarterGain = 1.;
+ StarterTorque = -1.;
+ StarterRPM = -1.;
// These are internal program variables
bBoostOverride = false;
bTakeoffBoost = false;
TakeoffBoost = 0.0; // Default to no extra takeoff-boost
+ BoostLossFactor = 0.0; // Default to free boost
+
int i;
for (i=0; i<FG_MAX_BOOST_SPEEDS; i++) {
RatedBoost[i] = 0.0;
// Read inputs from engine data file where present.
- if (el->FindElement("minmp")) // Should have ELSE statement telling default value used?
+ if (el->FindElement("minmp"))
MinManifoldPressure_inHg = el->FindElementValueAsNumberConvertTo("minmp","INHG");
if (el->FindElement("maxmp"))
MaxManifoldPressure_inHg = el->FindElementValueAsNumberConvertTo("maxmp","INHG");
+ if (el->FindElement("man-press-lag"))
+ ManifoldPressureLag = el->FindElementValueAsNumber("man-press-lag");
if (el->FindElement("displacement"))
Displacement = el->FindElementValueAsNumberConvertTo("displacement","IN3");
if (el->FindElement("maxhp"))
Ram_Air_Factor = el->FindElementValueAsNumber("ram-air-factor");
if (el->FindElement("cooling-factor"))
Cooling_Factor = el->FindElementValueAsNumber("cooling-factor");
+ if (el->FindElement("starter-rpm"))
+ StarterRPM = el->FindElementValueAsNumber("starter-rpm");
+ if (el->FindElement("starter-torque"))
+ StarterTorque = el->FindElementValueAsNumber("starter-torque");
if (el->FindElement("dynamic-fmep"))
FMEPDynamic= el->FindElementValueAsNumberConvertTo("dynamic-fmep","PA");
if (el->FindElement("static-fmep"))
BoostManual = (int)el->FindElementValueAsNumber("boostmanual");
if (el->FindElement("takeoffboost"))
TakeoffBoost = el->FindElementValueAsNumberConvertTo("takeoffboost", "PSI");
+ if (el->FindElement("boost-loss-factor"))
+ BoostLossFactor = el->FindElementValueAsNumber("boost-loss-factor");
if (el->FindElement("ratedboost1"))
RatedBoost[0] = el->FindElementValueAsNumberConvertTo("ratedboost1", "PSI");
if (el->FindElement("ratedboost2"))
}
}
- StarterHP = sqrt(MaxHP) * 0.4;
+
+ volumetric_efficiency_reduced = volumetric_efficiency;
+
+ if(StarterRPM < 0.) StarterRPM = 2*IdleRPM;
+ if(StarterTorque < 0)
+ StarterTorque = (MaxHP)*0.4; //just a wag.
+
displacement_SI = Displacement * in3tom3;
RatedMeanPistonSpeed_fps = ( MaxRPM * Stroke) / (360); // AKA 2 * (RPM/60) * ( Stroke / 12) or 2NS
base_property_name = CreateIndexedPropertyName("propulsion/engine", EngineNumber);
property_name = base_property_name + "/power-hp";
PropertyManager->Tie(property_name, &HP);
+ property_name = base_property_name + "/friction-hp";
+ PropertyManager->Tie(property_name, &StaticFriction_HP);
property_name = base_property_name + "/bsfc-lbs_hphr";
PropertyManager->Tie(property_name, &ISFC);
+ property_name = base_property_name + "/starter-norm";
+ PropertyManager->Tie(property_name, &StarterGain);
property_name = base_property_name + "/volumetric-efficiency";
PropertyManager->Tie(property_name, &volumetric_efficiency);
property_name = base_property_name + "/map-pa";
PropertyManager->Tie(property_name, this, &FGPiston::getOilPressure_psi);
property_name = base_property_name + "/egt-degF";
PropertyManager->Tie(property_name, this, &FGPiston::getExhaustGasTemp_degF);
+ if(BoostLossFactor > 0.0) {
+ property_name = base_property_name + "/boostloss-factor";
+ PropertyManager->Tie(property_name, &BoostLossFactor);
+ property_name = base_property_name + "/boostloss-hp";
+ PropertyManager->Tie(property_name, &BoostLossHP);
+ }
// Set up and sanity-check the turbo/supercharging configuration based on the input values.
if (TakeoffBoost > RatedBoost[0]) bTakeoffBoost = true;
Thruster->SetRPM(0.0);
RPM = 0.0;
OilPressure_psi = 0.0;
+ BoostLossHP = 0.;
}
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
int FGPiston::InitRunning(void)
{
Magnetos=3;
- in.MixtureCmd[EngineNumber] = in.PressureRatio/1.3;
- in.MixturePos[EngineNumber] = in.PressureRatio/1.3;
+ in.MixtureCmd[EngineNumber] = in.PressureRatio*1.3;
+ in.MixturePos[EngineNumber] = in.PressureRatio*1.3;
Thruster->SetRPM( 2.0*IdleRPM/Thruster->GetGearRatio() );
Running = true;
return 1;
if ((Magnetos == 1) || (Magnetos > 2)) Magneto_Left = true;
if (Magnetos > 1) Magneto_Right = true;
- // Assume we have fuel for now
- fuel = !Starved;
+// We will 'run' with any fuel flow. If there is not enough fuel to make power it will show in doEnginePower
+ fuel = FuelFlowRate > 0.0 ? 1 : 0;
// Check if we are turning the starter motor
if (Cranking != Starter) {
// This check saves .../cranking from getting updated every loop - they
// only update when changed.
Cranking = Starter;
- crank_counter = 0;
}
- if (Cranking) crank_counter++; //Check mode of engine operation
- if (!Running && spark && fuel) { // start the engine if revs high enough
- if (Cranking) {
- if ((RPM > IdleRPM*0.8) && (crank_counter > 175)) // Add a little delay to startup
- Running = true; // on the starter
- } else {
- if (RPM > IdleRPM*0.8) // This allows us to in-air start
- Running = true; // when windmilling
- }
- }
-
- // Cut the engine *power* - Note: the engine may continue to
- // spin if the prop is in a moving airstream
+ // Cut the engine *power* - Note: the engine will continue to
+ // spin depending on prop Ixx and freestream velocity
- if ( Running && (!spark || !fuel) ) Running = false;
-
- // Check for stalling (RPM = 0).
- if (Running) {
- if (RPM == 0) {
- Running = false;
- } else if ((RPM <= IdleRPM *0.8 ) && (Cranking)) {
- Running = false;
+ if ( Running ) {
+ if (!spark || !fuel) Running = false;
+ if (RPM < IdleRPM*0.8 ) Running = false;
+ } else { // !Running
+ if ( spark && fuel) { // start the engine if revs high enough
+ if (RPM > IdleRPM*0.8) // This allows us to in-air start
+ Running = true; // when windmilling
}
}
+
}
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
void FGPiston::doBoostControl(void)
{
- if(BoostManual) {
+ if(bBoostManual) {
if(BoostSpeed > BoostSpeeds-1) BoostSpeed = BoostSpeeds-1;
if(BoostSpeed < 0) BoostSpeed = 0;
} else {
* Inputs: p_amb, Throttle,
* MeanPistonSpeed_fps, dt
*
- * Outputs: MAP, ManifoldPressure_inHg, TMAP
+ * Outputs: MAP, ManifoldPressure_inHg, TMAP, BoostLossHP
*/
void FGPiston::doMAP(void)
double map_coefficient = Ze/(Ze+Z_airbox+Zt);
- // Add a one second lag to manifold pressure changes
- double dMAP=0;
- dMAP = (TMAP - p_ram * map_coefficient) * TotalDeltaT;
+ // Add a variable lag to manifold pressure changes
+ double dMAP=(TMAP - p_ram * map_coefficient);
+ if (ManifoldPressureLag > TotalDeltaT) dMAP *= TotalDeltaT/ManifoldPressureLag;
TMAP -=dMAP;
double boost_factor = (( BoostMul[BoostSpeed] - 1 ) / RatedRPM[BoostSpeed] ) * RPM + 1;
MAP = TMAP * boost_factor;
// Now clip the manifold pressure to BCV or Wastegate setting.
- if (bTakeoffPos) {
- if (MAP > TakeoffMAP[BoostSpeed]) MAP = TakeoffMAP[BoostSpeed];
- } else {
- if (MAP > RatedMAP[BoostSpeed]) MAP = RatedMAP[BoostSpeed];
+ if(!bBoostOverride) {
+ if (bTakeoffPos) {
+ if (MAP > TakeoffMAP[BoostSpeed]) MAP = TakeoffMAP[BoostSpeed];
+ } else {
+ if (MAP > RatedMAP[BoostSpeed]) MAP = RatedMAP[BoostSpeed];
+ }
}
} else {
MAP = TMAP;
}
+ if( BoostLossFactor > 0.0 )
+ {
+ double gamma = 1.414; // specific heat constants
+ double Nstage = 1; // Nstage is the number of boost stages.
+ BoostLossHP = ((Nstage * TMAP * v_dot_air * gamma) / (gamma - 1)) * (pow((MAP/TMAP),((gamma-1)/(Nstage * gamma))) - 1) * BoostLossFactor / 745.7 ; // 745.7 convert watt to hp
+ } else {
+ BoostLossHP = 0;
+ }
+
// And set the value in American units as well
ManifoldPressure_inHg = MAP / inhgtopa;
}
*
* TODO: Model inlet manifold air temperature.
*
- * Outputs: rho_air, m_dot_air
+ * Outputs: rho_air, m_dot_air, volumetric_efficiency_reduced
*/
void FGPiston::doAirFlow(void)
double gamma = 1.3; // specific heat constants
// loss of volumentric efficiency due to difference between MAP and exhaust pressure
// Eq 6-10 from The Internal Combustion Engine - Charles Taylor Vol 1
- double ve =((gamma-1)/gamma) +( CompressionRatio -(p_amb/MAP))/(gamma*( CompressionRatio - 1));
+ double mratio = MAP < 1. ? CompressionRatio : p_amb/MAP;
+ if (mratio > CompressionRatio) mratio = CompressionRatio;
+ double ve =((gamma-1)/gamma) +( CompressionRatio -(mratio))/(gamma*( CompressionRatio - 1));
rho_air = p_amb / (R_air * T_amb);
double swept_volume = (displacement_SI * (RPM/60)) / 2;
- double v_dot_air = swept_volume * volumetric_efficiency *ve;
+ volumetric_efficiency_reduced = volumetric_efficiency *ve;
+ v_dot_air = swept_volume * volumetric_efficiency_reduced;
double rho_air_manifold = MAP / (R_air * T_amb);
m_dot_air = v_dot_air * rho_air_manifold;
void FGPiston::doEnginePower(void)
{
- IndicatedHorsePower = 0;
+ IndicatedHorsePower = -StaticFriction_HP;
FMEP = 0;
if (Running) {
- // FIXME: this needs to be generalized
- double ME, percent_RPM, power; // Convienience term for use in the calculations
+ double ME, power; // Convienience term for use in the calculations
ME = Mixture_Efficiency_Correlation->GetValue(m_dot_fuel/m_dot_air);
- percent_RPM = RPM/MaxRPM;
// Guestimate engine friction losses from Figure 4.4 of "Engines: An Introduction", John Lumley
FMEP = (-FMEPDynamic * MeanPistonSpeed_fps * fttom - FMEPStatic);
if ( Magnetos != 3 ) power *= SparkFailDrop;
- IndicatedHorsePower = (FuelFlow_pph / ISFC )* ME * power;
+ IndicatedHorsePower = (FuelFlow_pph / ISFC )* ME * power - StaticFriction_HP; //FIXME static friction should depend on oil temp and configuration;
} else {
// Power output when the engine is not running
+ double torque, k_torque, rpm; // Convienience term for use in the calculations
+
+ rpm = RPM < 1.0 ? 1.0 : RPM;
if (Cranking) {
- if (RPM < 10) {
- IndicatedHorsePower = StarterHP;
- } else if (RPM < IdleRPM*0.8) {
- IndicatedHorsePower = StarterHP + ((IdleRPM*0.8 - RPM) / 8.0);
- // This is a guess - would be nice to find a proper starter moter torque curve
- } else {
- IndicatedHorsePower = StarterHP;
- }
- }
+ if(RPM<StarterRPM) k_torque = 1.0-RPM/(StarterRPM);
+ else k_torque = 0;
+ torque = StarterTorque*k_torque*StarterGain;
+ IndicatedHorsePower = torque * rpm / 5252;
+ }
}
// Constant is (1/2) * 60 * 745.7
// (1/2) convert cycles, 60 minutes to seconds, 745.7 watts to hp.
double pumping_hp = ((PMEP + FMEP) * displacement_SI * RPM)/(Cycles*22371);
- HP = IndicatedHorsePower + pumping_hp - StaticFriction_HP; //FIXME static friction should depend on oil temp and configuration
+HP = IndicatedHorsePower + pumping_hp - BoostLossHP;
// cout << "pumping_hp " <<pumping_hp << FMEP << PMEP <<endl;
PctPower = HP / MaxHP ;
// cout << "Power = " << HP << " RPM = " << RPM << " Running = " << Running << " Cranking = " << Cranking << endl;
cout << " Bore: " << Bore << endl;
cout << " Stroke: " << Stroke << endl;
cout << " Cylinders: " << Cylinders << endl;
- cout << " Cylinders Head Mass: " <<CylinderHeadMass << endl;
+ cout << " Cylinders Head Mass: " << CylinderHeadMass << endl;
cout << " Compression Ratio: " << CompressionRatio << endl;
cout << " MaxHP: " << MaxHP << endl;
cout << " Cycles: " << Cycles << endl;
cout << " Dynamic FMEP Factor: " << FMEPDynamic << endl;
cout << " Static FMEP Factor: " << FMEPStatic << endl;
+ cout << " Starter Motor Torque: " << StarterTorque << endl;
+ cout << " Starter Motor RPM: " << StarterRPM << endl;
+
cout << endl;
cout << " Combustion Efficiency table:" << endl;
Lookup_Combustion_Efficiency->Print();
DEFINITIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#define ID_PISTON "$Id: FGPiston.h,v 1.32 2011/10/11 16:16:16 jentron Exp $";
+#define ID_PISTON "$Id: FGPiston.h,v 1.35 2012/04/07 01:50:54 jentron Exp $";
#define FG_MAX_BOOST_SPEEDS 3
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
<air-intake-impedance-factor> {number} </air-intake-impedance-factor>
<ram-air-factor> {number} </ram-air-factor>
<cooling-factor> {number} </cooling-factor>
+ <man-press-lag> {number} </man-press-lag>
+ <starter-torque> {number} </starter-torque>
+ <starter-rpm> {number} </starter-rpm>
<cylinder-head-mass unit="{KG | LBS}"> {number} </cylinder-head-mass>
<bsfc unit="{LBS/HP*HR | "KG/KW*HR"}"> {number} </bsfc>
<volumetric-efficiency> {number} </volumetric-efficiency>
<numboostspeeds> {number} </numboostspeeds>
<boostoverride> {0 | 1} </boostoverride>
<boostmanual> {0 | 1} </boostmanual>
+ <boost-loss-factor> {number} </boost-loss-factor>
<ratedboost1 unit="{INHG | PA | ATM}"> {number} </ratedboost1>
<ratedpower1 unit="{HP | WATTS}"> {number} </ratedpower1>
<ratedrpm1> {number} </ratedrpm1>
- \b maxmp - this value is the nominal maximum manifold pressure at sea-level
without boost. Along with maxrpm it determines the resistance of the
aircraft intake system. Overridden by air-intake-impedance-factor
+- \b man-press-lag - Delay in seconds for manifold pressure changes to take effect
+- \b starter-torque - A value specifing the zero RPM torque in lb*ft the starter motor
+ provides. Current default value is 40% of the horse power value.
+- \b starter-rpm - A value specifing the maximum RPM the unloaded starter motor
+ can achieve. Loads placed on the engine by the propeller and throttle will
+ further limit RPM achieved in practice.
- \b idlerpm - this value affects the throttle fall off and the engine stops
running if it is slowed below 80% of this value. The engine starts
running when it reaches 80% of this value.
supercharger speeds. Merlin XII had 1 speed, Merlin 61 had 2, a late
Griffon engine apparently had 3. No known engine more than 3, although
some German engines had continuously variable-speed superchargers.
-- \b boostoverride - unused
+- \b boostoverride - whether or not to clip output to the wastegate value
+- \b boost-loss-factor - zero (or not present) for 'free' supercharging. A value entered
+ will be used as a multiplier to the power required to compress the input air. Typical
+ value should be 1.15 to 1.20.
- \b boostmanual - whether a multispeed supercharger will manually or
automatically shift boost speeds. On manual shifting the boost speeds is
accomplished by controlling the property propulsion/engine/boostspeed.
@author David Megginson (initial porting and additional code)
@author Ron Jensen (additional engine code)
@see Taylor, Charles Fayette, "The Internal Combustion Engine in Theory and Practice"
- @version $Id: FGPiston.h,v 1.32 2011/10/11 16:16:16 jentron Exp $
+ @version $Id: FGPiston.h,v 1.35 2012/04/07 01:50:54 jentron Exp $
*/
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
int crank_counter;
double IndicatedHorsePower;
+ double IndicatedPower;
double PMEP;
double FMEP;
double FMEPDynamic;
double FMEPStatic;
+ double T_Intake;
void doEngineStartup(void);
void doBoostControl(void);
double MinManifoldPressure_inHg; // Inches Hg
double MaxManifoldPressure_inHg; // Inches Hg
double MaxManifoldPressure_Percent; // MaxManifoldPressure / 29.92
+ double ManifoldPressureLag; // Manifold Pressure delay in seconds.
double Displacement; // cubic inches
double displacement_SI; // cubic meters
double MaxHP; // horsepower
double RatedMeanPistonSpeed_fps; // ft/sec derived from MaxRPM and stroke.
double Ram_Air_Factor; // number
- double StarterHP; // initial horsepower of starter motor
+ double StarterTorque;// Peak Torque of the starter motor
+ double StarterRPM; // Peak RPM of the starter motor
+ double StarterGain; // control the torque of the starter motor.
int BoostSpeeds; // Number of super/turbocharger boost speeds - zero implies no turbo/supercharging.
int BoostSpeed; // The current boost-speed (zero-based).
bool Boosted; // Set true for boosted engine.
double RatedMAP[FG_MAX_BOOST_SPEEDS]; // Rated manifold absolute pressure [Pa] (BCV clamp)
double TakeoffMAP[FG_MAX_BOOST_SPEEDS]; // Takeoff setting manifold absolute pressure [Pa] (BCV clamp)
double BoostSwitchHysteresis; // Pa.
+ double BoostLossFactor; // multiplier for HP consumed by the supercharger
double minMAP; // Pa
double maxMAP; // Pa
//
double rho_air;
double volumetric_efficiency;
+ double volumetric_efficiency_reduced;
double map_coefficient;
double m_dot_air;
+ double v_dot_air;
double equivalence_ratio;
double m_dot_fuel;
double HP;
+ double BoostLossHP;
double combustion_efficiency;
double ExhaustGasTemp_degK;
double EGT_degC;
namespace JSBSim {
-static const char *IdSrc = "$Id: FGPropeller.cpp,v 1.41 2011/11/17 21:07:30 jentron Exp $";
+static const char *IdSrc = "$Id: FGPropeller.cpp,v 1.44 2012/04/29 13:10:46 bcoconni Exp $";
static const char *IdHdr = ID_PROPELLER;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
// Induced velocity in the propeller disk area. This formula is obtained
// from momentum theory - see B. W. McCormick, "Aerodynamics, Aeronautics,
// and Flight Mechanics" 1st edition, eqn. 6.15 (propeller analysis chapter).
- // Vinduced = 0.5 * (-Vel + sqrt(Vel*Vel + 2.0*Thrust/(rho*Area)))
// Since Thrust and Vel can both be negative we need to adjust this formula
// To handle sign (direction) separately from magnitude.
double Vel2sum = Vel*abs(Vel) + 2.0*Thrust/(rho*Area);
-
+
if( Vel2sum > 0.0)
Vinduced = 0.5 * (-Vel + sqrt(Vel2sum));
else
Vinduced = 0.5 * (-Vel - sqrt(-Vel2sum));
+ // We need to drop the case where the downstream velocity is opposite in
+ // direction to the aircraft velocity. For example, in such a case, the
+ // direction of the airflow on the tail would be opposite to the airflow on
+ // the wing tips. When such complicated airflows occur, the momentum theory
+ // breaks down and the formulas above are no longer applicable
+ // (see H. Glauert, "The Elements of Airfoil and Airscrew Theory",
+ // 2nd edition, §16.3, pp. 219-221)
+
+ if ((Vel+2.0*Vinduced)*Vel < 0.0)
+ Vinduced = 0.0; // We cannot calculate the induced velocity so let's assume it is zero.
+
// P-factor is simulated by a shift of the acting location of the thrust.
// The shift is a multiple of the angle between the propeller shaft axis
// and the relative wind that goes through the propeller disk.
namespace JSBSim {
-static const char *IdSrc = "$Id: FGRocket.cpp,v 1.26 2011/08/04 13:45:42 jberndt Exp $";
+static const char *IdSrc = "$Id: FGRocket.cpp,v 1.27 2012/04/08 15:19:08 jberndt Exp $";
static const char *IdHdr = ID_ROCKET;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
FGRocket::FGRocket(FGFDMExec* exec, Element *el, int engine_number, struct Inputs& input)
- : FGEngine(exec, el, engine_number, input)
+ : FGEngine(exec, el, engine_number, input), isp_function(0L)
{
Type = etRocket;
Element* thrust_table_element = 0;
TotalPropellantExpended = 0.0;
FuelFlowRate = FuelExpended = 0.0;
OxidizerFlowRate = OxidizerExpended = 0.0;
- SLOxiFlowMax = SLFuelFlowMax = 0.0;
+ SLOxiFlowMax = SLFuelFlowMax = PropFlowMax = 0.0;
+ MxR = 0.0;
BuildupTime = 0.0;
It = 0.0;
ThrustVariation = 0.0;
MinThrottle = 0.0;
MaxThrottle = 1.0;
- if (el->FindElement("isp"))
+ string base_property_name = CreateIndexedPropertyName("propulsion/engine", EngineNumber);
+
+ std::stringstream strEngineNumber;
+ strEngineNumber << EngineNumber;
+
+ Element* isp_el = el->FindElement("isp");
+ Element* isp_func_el=0;
+
+ bindmodel(); // Bind model properties first, since they might be needed in functions.
+
+ // Specific impulse may be specified as a constant value or as a function - perhaps as a function of mixture ratio.
+ if (isp_el) {
+ isp_func_el = isp_el->FindElement("function");
+ if (isp_func_el) {
+ isp_function = new FGFunction(exec->GetPropertyManager(),isp_func_el, strEngineNumber.str());
+ } else {
Isp = el->FindElementValueAsNumber("isp");
+ }
+ } else {
+ throw("Specific Impulse <isp> must be specified for a rocket engine");
+ }
+
if (el->FindElement("builduptime"))
BuildupTime = el->FindElementValueAsNumber("builduptime");
if (el->FindElement("maxthrottle"))
MaxThrottle = el->FindElementValueAsNumber("maxthrottle");
if (el->FindElement("minthrottle"))
MinThrottle = el->FindElementValueAsNumber("minthrottle");
- if (el->FindElement("slfuelflowmax"))
+
+ if (el->FindElement("slfuelflowmax")) {
SLFuelFlowMax = el->FindElementValueAsNumberConvertTo("slfuelflowmax", "LBS/SEC");
- if (el->FindElement("sloxiflowmax"))
+ if (el->FindElement("sloxiflowmax")) {
SLOxiFlowMax = el->FindElementValueAsNumberConvertTo("sloxiflowmax", "LBS/SEC");
+ }
+ PropFlowMax = SLOxiFlowMax + SLFuelFlowMax;
+ MxR = SLOxiFlowMax/SLFuelFlowMax;
+ } else if (el->FindElement("propflowmax")) {
+ PropFlowMax = el->FindElementValueAsNumberConvertTo("propflowmax", "LBS/SEC");
+ // Mixture ratio may be specified here, but it can also be specified as a function or via property
+ if (el->FindElement("mixtureratio")) {
+ MxR = el->FindElementValueAsNumber("mixtureratio");
+ }
+ }
+ if (isp_function) Isp = isp_function->GetValue(); // cause Isp function to be executed if present.
// If there is a thrust table element, this is a solid propellant engine.
thrust_table_element = el->FindElement("thrust_table");
if (thrust_table_element) {
}
}
- bindmodel();
Debug(0);
}
PropellantFlowRate = (FuelExpended + OxidizerExpended)/in.TotalDeltaT;
TotalPropellantExpended += FuelExpended + OxidizerExpended;
+ // If Isp has been specified as a function, override the value of Isp to that, otherwise
+ // assume a constant value is given.
+ if (isp_function) Isp = isp_function->GetValue();
// If there is a thrust table, it is a function of propellant burned. The
// engine is started when the throttle is advanced to 1.0. After that, it
FuelFlowRate = VacThrust/Isp; // This calculates wdot (weight flow rate in lbs/sec)
FuelFlowRate /= (1 + TotalIspVariation);
} else {
- FuelFlowRate = SLFuelFlowMax*PctPower;
+ SLFuelFlowMax = PropFlowMax / (1 + MxR);
+ FuelFlowRate = SLFuelFlowMax * PctPower;
}
FuelExpended = FuelFlowRate * in.TotalDeltaT; // For this time step ...
double FGRocket::CalcOxidizerNeed(void)
{
+ SLOxiFlowMax = PropFlowMax * MxR / (1 + MxR);
OxidizerFlowRate = SLOxiFlowMax * PctPower;
OxidizerExpended = OxidizerFlowRate * in.TotalDeltaT;
return OxidizerExpended;
} else { // Liquid rocket motor
property_name = base_property_name + "/oxi-flow-rate-pps";
PropertyManager->Tie( property_name.c_str(), this, &FGRocket::GetOxiFlowRate);
+ property_name = base_property_name + "/mixture-ratio";
+ PropertyManager->Tie( property_name.c_str(), this, &FGRocket::GetMixtureRatio,
+ &FGRocket::SetMixtureRatio);
+ property_name = base_property_name + "/isp";
+ PropertyManager->Tie( property_name.c_str(), this, &FGRocket::GetIsp,
+ &FGRocket::SetIsp);
}
}
#include "FGEngine.h"
#include "math/FGTable.h"
+#include "math/FGFunction.h"
#include "input_output/FGXMLElement.h"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
DEFINITIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#define ID_ROCKET "$Id: FGRocket.h,v 1.17 2011/08/04 13:45:42 jberndt Exp $"
+#define ID_ROCKET "$Id: FGRocket.h,v 1.18 2012/04/08 15:19:08 jberndt Exp $"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FORWARD DECLARATIONS
fuel begins burning and thrust is provided.
@author Jon S. Berndt
- $Id: FGRocket.h,v 1.17 2011/08/04 13:45:42 jberndt Exp $
+ $Id: FGRocket.h,v 1.18 2012/04/08 15:19:08 jberndt Exp $
@see FGNozzle,
FGThruster,
FGForce,
double GetOxiFlowRate(void) const {return OxidizerFlowRate;}
+ double GetMixtureRatio(void) const {return MxR;}
+
+ double GetIsp(void) const {return Isp;}
+
+ void SetMixtureRatio(double mix) {MxR = mix;}
+
+ void SetIsp(double isp) {Isp = isp;}
+
std::string GetEngineLabels(const std::string& delimiter);
std::string GetEngineValues(const std::string& delimiter);
double OxidizerExpended;
double TotalPropellantExpended;
double SLOxiFlowMax;
+ double PropFlowMax;
double OxidizerFlowRate;
double PropellantFlowRate;
bool Flameout;
double BuildupTime;
FGTable* ThrustTable;
+ FGFunction* isp_function;
void Debug(int from);
};
01/10/11 T.Kreitler changed to single rotor model
03/06/11 T.Kreitler added brake, clutch, and experimental free-wheeling-unit,
reasonable estimate for inflowlag
+02/05/12 T.Kreitler brake, clutch, and FWU now in FGTransmission,
+ downwash angles relate to shaft orientation
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
INCLUDES
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#include <iostream>
#include <sstream>
#include "FGRotor.h"
-#include "input_output/FGXMLElement.h"
#include "models/FGMassBalance.h"
#include "models/FGPropulsion.h" // to get the GearRatio from a linked rotor
using std::cerr;
+using std::cout;
using std::endl;
using std::ostringstream;
-using std::cout;
namespace JSBSim {
-static const char *IdSrc = "$Id: FGRotor.cpp,v 1.18 2011/10/15 21:30:28 jentron Exp $";
+static const char *IdSrc = "$Id: FGRotor.cpp,v 1.20 2012/03/18 15:48:36 jentron Exp $";
static const char *IdHdr = ID_ROTOR;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
CLASS IMPLEMENTATION
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-
-//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-
-// Note: The FGTransmission class is currently carried 'pick-a-pack' by
-// FGRotor.
-
-FGTransmission::FGTransmission(FGFDMExec *exec, int num) :
- FreeWheelTransmission(1.0),
- ThrusterMoment(1.0), EngineMoment(1.0), EngineFriction(0.0),
- ClutchCtrlNorm(1.0), BrakeCtrlNorm(0.0), MaxBrakePower(0.0),
- EngineRPM(0.0), ThrusterRPM(0.0)
-{
- double dt;
- PropertyManager = exec->GetPropertyManager();
- dt = exec->GetDeltaT();
-
- // avoid too abrupt changes in transmission
- FreeWheelLag = Filter(200.0,dt);
- BindModel(num);
-}
-
-//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-
-FGTransmission::~FGTransmission(){
-}
-
-//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-
-// basically P = Q*w and Q_Engine + (-Q_Rotor) = J * dw/dt, J = Moment
-//
-void FGTransmission::Calculate(double EnginePower, double ThrusterTorque, double dt) {
-
- double coupling = 1.0, coupling_sq = 1.0;
- double fw_mult = 1.0;
-
- double d_omega = 0.0, engine_d_omega = 0.0, thruster_d_omega = 0.0; // relative changes
-
- double engine_omega = rpm_to_omega(EngineRPM);
- double safe_engine_omega = engine_omega < 1e-1 ? 1e-1 : engine_omega;
- double engine_torque = EnginePower / safe_engine_omega;
-
- double thruster_omega = rpm_to_omega(ThrusterRPM);
- double safe_thruster_omega = thruster_omega < 1e-1 ? 1e-1 : thruster_omega;
-
- engine_torque -= EngineFriction / safe_engine_omega;
- ThrusterTorque += Constrain(0.0, BrakeCtrlNorm, 1.0) * MaxBrakePower / safe_thruster_omega;
-
- // would the FWU release ?
- engine_d_omega = engine_torque/EngineMoment * dt;
- thruster_d_omega = - ThrusterTorque/ThrusterMoment * dt;
-
- if ( thruster_omega+thruster_d_omega > engine_omega+engine_d_omega ) {
- // don't drive the engine
- FreeWheelTransmission = 0.0;
- } else {
- FreeWheelTransmission = 1.0;
- }
-
- fw_mult = FreeWheelLag.execute(FreeWheelTransmission);
- coupling = fw_mult * Constrain(0.0, ClutchCtrlNorm, 1.0);
-
- if (coupling < 0.999999) { // are the separate calculations needed ?
- // assume linear transfer
- engine_d_omega =
- (engine_torque - ThrusterTorque*coupling)/(ThrusterMoment*coupling + EngineMoment) * dt;
- thruster_d_omega =
- (engine_torque*coupling - ThrusterTorque)/(ThrusterMoment + EngineMoment*coupling) * dt;
-
- EngineRPM += omega_to_rpm(engine_d_omega);
- ThrusterRPM += omega_to_rpm(thruster_d_omega);
-
- // simulate friction in clutch
- coupling_sq = coupling*coupling;
- EngineRPM = (1.0-coupling_sq) * EngineRPM + coupling_sq * 0.02 * (49.0*EngineRPM + ThrusterRPM);
- ThrusterRPM = (1.0-coupling_sq) * ThrusterRPM + coupling_sq * 0.02 * (EngineRPM + 49.0*ThrusterRPM);
-
- // enforce equal rpm
- if ( fabs(EngineRPM-ThrusterRPM) < 1e-3 ) {
- EngineRPM = ThrusterRPM = 0.5 * (EngineRPM+ThrusterRPM);
- }
- } else {
- d_omega = (engine_torque - ThrusterTorque)/(ThrusterMoment + EngineMoment) * dt;
- EngineRPM = ThrusterRPM += omega_to_rpm(d_omega);
- }
-
- // nothing will turn backward
- if (EngineRPM < 0.0 ) EngineRPM = 0.0;
- if (ThrusterRPM < 0.0 ) ThrusterRPM = 0.0;
-
-}
-
-//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-
-bool FGTransmission::BindModel(int num)
-{
- string property_name, base_property_name;
- base_property_name = CreateIndexedPropertyName("propulsion/engine", num);
-
- property_name = base_property_name + "/brake-ctrl-norm";
- PropertyManager->Tie( property_name.c_str(), this, &FGTransmission::GetBrakeCtrl, &FGTransmission::SetBrakeCtrl);
- property_name = base_property_name + "/free-wheel-transmission";
- PropertyManager->Tie( property_name.c_str(), this, &FGTransmission::GetFreeWheelTransmission);
-
- return true;
-}
-
-
-
-//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-
// Constructor
FGRotor::FGRotor(FGFDMExec *exec, Element* rotor_element, int num)
BladeChord(0.0), LiftCurveSlope(0.0), BladeTwist(0.0), HingeOffset(0.0),
BladeFlappingMoment(0.0), BladeMassMoment(0.0), PolarMoment(0.0),
InflowLag(0.0), TipLossB(0.0),
- GroundEffectExp(0.0), GroundEffectShift(0.0),
+ GroundEffectExp(0.0), GroundEffectShift(0.0), GroundEffectScaleNorm(1.0),
LockNumberByRho(0.0), Solidity(0.0), // derived parameters
RPM(0.0), Omega(0.0), // dynamic values
beta_orient(0.0),
{
FGColumnVector3 location(0.0, 0.0, 0.0), orientation(0.0, 0.0, 0.0);
Element *thruster_element;
+ double engine_power_est = 0.0;
// initialise/set remaining variables
SetTransformType(FGForce::tCustom);
SetLocation(location);
SetAnglesToBody(orientation);
- InvTransform = Transform().Transposed();
+ InvTransform = Transform().Transposed(); // body to custom/native
// wire controls
ControlMap = eMainCtrl;
}
}
- // configure the rotor parameters
- Configure(rotor_element);
+ // process rotor parameters
+ engine_power_est = Configure(rotor_element);
- // configure the transmission parameters
+ // setup transmission if needed
if (!ExternalRPM) {
- Transmission = new FGTransmission(exec, num);
- Transmission->SetMaxBrakePower(MaxBrakePower);
-
- if (rotor_element->FindElement("gearloss")) {
- GearLoss = rotor_element->FindElementValueAsNumberConvertTo("gearloss","HP");
- GearLoss *= hptoftlbssec;
- }
+ Transmission = new FGTransmission(exec, num, dt);
- if (GearLoss<1e-6) { // TODO, allow 0 ?
- double ehp = 0.5 * BladeNum*BladeChord*Radius*Radius; // guess engine power
- //cout << "# guessed engine power: " << ehp << endl;
- GearLoss = 0.0025 * ehp * hptoftlbssec;
- }
- Transmission->SetEngineFriction(GearLoss);
+ Transmission->SetThrusterMoment(PolarMoment);
// The MOI sensed behind the gear ( MOI_engine*sqr(GearRatio) ).
- if (rotor_element->FindElement("gearmoment")) {
- GearMoment = rotor_element->FindElementValueAsNumberConvertTo("gearmoment","SLUG*FT2");
- }
-
- if (GearMoment<1e-2) { // TODO, need better check for lower limit
- GearMoment = 0.1*PolarMoment;
- }
+ GearMoment = ConfigValueConv(rotor_element, "gearmoment", 0.1*PolarMoment, "SLUG*FT2");
+ GearMoment = Constrain(1e-6, GearMoment, 1e9);
Transmission->SetEngineMoment(GearMoment);
- Transmission->SetThrusterMoment(PolarMoment);
+ Transmission->SetMaxBrakePower(MaxBrakePower);
+
+ GearLoss = ConfigValueConv(rotor_element, "gearloss", 0.0025 * engine_power_est, "HP");
+ GearLoss = Constrain(0.0, GearLoss, 1e9);
+ GearLoss *= hptoftlbssec;
+ Transmission->SetEngineFriction(GearLoss);
+
}
- // shaft representation - a rather simple transform,
+ // shaft representation - a rather simple transform,
// but using a matrix is safer.
TboToHsr.InitMatrix( 0.0, 0.0, 1.0,
0.0, 1.0, 0.0,
Debug(1);
}
-
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-// 5in1: value-fetch-convert-default-return function
+// 5in1: value-fetch-convert-default-return function
-double FGRotor::ConfigValueConv( Element* el, const string& ename, double default_val,
+double FGRotor::ConfigValueConv( Element* el, const string& ename, double default_val,
const string& unit, bool tell)
{
// 1. read configuration and try to fill holes, ymmv
// 2. calculate derived parameters
-void FGRotor::Configure(Element* rotor_element)
+double FGRotor::Configure(Element* rotor_element)
{
- double estimate;
+ double estimate, engine_power_est=0.0;
const bool yell = true;
const bool silent = false;
- Radius = 0.5 * ConfigValueConv(rotor_element, "diameter", 42.0, "FT", yell);
+ Radius = 0.5 * ConfigValueConv(rotor_element, "diameter", 42.0, "FT", yell);
Radius = Constrain(1e-3, Radius, 1e9);
BladeNum = (int) ConfigValue(rotor_element, "numblades", 3 , yell);
GearRatio = ConfigValue(rotor_element, "gearratio", 1.0, yell);
+ GearRatio = Constrain(1e-9, GearRatio, 1e9);
// make sure that v_tip (omega*r) is below 0.7mach ~ 750ft/s
estimate = (750.0/Radius)/(2.0*M_PI) * 60.0; // 7160/Radius
estimate = sqr(BladeChord) * sqr(Radius - HingeOffset) * 0.57;
BladeFlappingMoment = ConfigValueConv(rotor_element, "flappingmoment", estimate, "SLUG*FT2");
- BladeFlappingMoment = Constrain(1.0e-6, BladeFlappingMoment, 1e9);
+ BladeFlappingMoment = Constrain(1e-9, BladeFlappingMoment, 1e9);
// guess mass from moment of a thin stick, and multiply by the blades cg distance
estimate = ( 3.0 * BladeFlappingMoment / sqr(Radius) ) * (0.45 * Radius) ;
BladeMassMoment = ConfigValue(rotor_element, "massmoment", estimate); // unit is slug-ft
- BladeMassMoment = Constrain(0.001, BladeMassMoment, 1e9);
+ BladeMassMoment = Constrain(1e-9, BladeMassMoment, 1e9);
estimate = 1.1 * BladeFlappingMoment * BladeNum;
PolarMoment = ConfigValueConv(rotor_element, "polarmoment", estimate, "SLUG*FT2");
- PolarMoment = Constrain(1e-6, PolarMoment, 1e9);
+ PolarMoment = Constrain(1e-9, PolarMoment, 1e9);
// "inflowlag" is treated further down.
TipLossB = ConfigValue(rotor_element, "tiplossfactor", 1.0, silent);
- estimate = 0.01 * PolarMoment ; // guesses for huey, bo105 20-30hp
+ // estimate engine power (bit of a pity, cause our caller already knows)
+ engine_power_est = 0.5 * BladeNum*BladeChord*Radius*Radius;
+
+ estimate = engine_power_est / 30.0;
MaxBrakePower = ConfigValueConv(rotor_element, "maxbrakepower", estimate, "HP");
MaxBrakePower *= hptoftlbssec;
- // ground effect
- if (rotor_element->FindElement("cgroundeffect")) {
- double cge,gee;
- cge = rotor_element->FindElementValueAsNumber("cgroundeffect");
- cge = Constrain(1e-9, cge, 1.0);
- gee = 1.0 / ( 2.0*Radius * cge );
- cerr << "# *** 'cgroundeffect' is defunct." << endl;
- cerr << "# *** use 'groundeffectexp' with: " << gee << endl;
- }
-
GroundEffectExp = ConfigValue(rotor_element, "groundeffectexp", 0.0);
GroundEffectShift = ConfigValueConv(rotor_element, "groundeffectshift", 0.0, "FT");
estimate = 16.0/(LockNumberByRho*rho * omega_tmp ); // 16/(gamma*Omega)
// printf("# Est. InflowLag: %f\n", estimate);
InflowLag = ConfigValue(rotor_element, "inflowlag", estimate, yell);
- InflowLag = Constrain(1.0e-6, InflowLag, 2.0);
+ InflowLag = Constrain(1e-6, InflowLag, 2.0);
- return;
+ return engine_power_est;
} // Configure
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FGColumnVector3 FGRotor::fus_angvel_body2ca( const FGColumnVector3 &pqr)
{
- FGColumnVector3 av_s_fus, av_w_fus;
+ FGColumnVector3 av_s_fus, av_w_fus;
// for comparison:
- // av_s_fus = BodyToShaft * pqr; /SH79/
+ // av_s_fus = BodyToShaft * pqr; /SH79/
// BodyToShaft = TboToHsr * InvTransform
av_s_fus = TboToHsr * InvTransform * pqr;
// reduction of inflow if the helicopter is close to the ground, yielding to
// higher thrust, see /TA77/ eqn(10a).
-void FGRotor::calc_flow_and_thrust( double theta_0, double Uw, double Ww,
+void FGRotor::calc_flow_and_thrust( double theta_0, double Uw, double Ww,
double flow_scale)
{
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-// The coning angle doesn't apply for teetering rotors, but calculating
-// doesn't hurt. /SH79/ eqn(29)
+// Two blade teetering rotors are often 'preconed' to a fixed angle, but the
+// calculated value is pretty close to the real one. /SH79/ eqn(29)
void FGRotor::calc_coning_angle(double theta_0)
{
// estimate blade drag
double delta_dr = 0.009 + 0.3*sqr(6.0*C_T/(LiftCurveSlope*Solidity));
- Torque = rho * BladeNum * BladeChord * delta_dr * sqr(Omega*Radius) * R[2] *
+ Torque = rho * BladeNum * BladeChord * delta_dr * sqr(Omega*Radius) * R[2] *
(1.0+4.5*sqr(mu))/8.0
- (Thrust*lambda + H_drag*mu)*Radius;
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+// Get the downwash angles with respect to the shaft axis.
+// Given a 'regular' main rotor, the angles are zero when the downwash points
+// down, positive theta values mean that the downwash turns towards the nose,
+// and positive phi values mean it turns to the left side. (Note: only airspeed
+// is transformed, the rotational speed contribution is ignored.)
+
+void FGRotor::calc_downwash_angles()
+{
+ FGColumnVector3 v_shaft;
+ v_shaft = TboToHsr * InvTransform * in.AeroUVW;
+
+ theta_downwash = atan2( -v_shaft(eU), v_induced - v_shaft(eW)) + a1s;
+ phi_downwash = atan2( v_shaft(eV), v_induced - v_shaft(eW)) + b1s;
+
+ return;
+}
+
+//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
// transform rotor forces from control axes to shaft axes, and express
// in body axes /SH79/ eqn(40,41)
FGColumnVector3 F_s(
- H_drag*cos(beta_orient) - J_side*sin(beta_orient) + Thrust*b_ic,
- H_drag*sin(beta_orient) + J_side*cos(beta_orient) + Thrust*a_ic,
- - Thrust);
+ - Thrust);
return HsrToTbo * F_s;
}
// fetch needed values from environment
rho = in.Density; // slugs/ft^3.
double h_agl_ft = in.H_agl;
+
// update InvTransform, the rotor orientation could have been altered
InvTransform = Transform().Transposed();
B_IC = LongitudinalCtrl;
theta_col = CollectiveCtrl;
- // ground effect
+ // optional ground effect, a ge_factor of 1.0 gives no effect
+ // and 0.5 yields to maximal influence.
if (GroundEffectExp > 1e-5) {
if (h_agl_ft<0.0) h_agl_ft = 0.0; // clamp
filtered_hagl = damp_hagl.execute(h_agl_ft) + GroundEffectShift;
// actual/nominal factor avoids absurd scales at startup
- ge_factor -= exp(-filtered_hagl*GroundEffectExp) * (RPM / NominalRPM);
- if (ge_factor<0.5) ge_factor=0.5; // clamp
+ ge_factor -= GroundEffectScaleNorm *
+ ( exp(-filtered_hagl*GroundEffectExp) * (RPM / NominalRPM) );
+ ge_factor = Constrain(0.5, ge_factor, 1.0);
}
- // all set, start calculations
+ // all set, start calculations ...
vHub_ca = hub_vel_body2ca(in.AeroUVW, in.AeroPQR, A_IC, B_IC);
calc_torque(theta_col);
- // Fixme: only valid for a 'decent' rotor
- theta_downwash = atan2( -in.AeroUVW(eU), v_induced - in.AeroUVW(eW));
- phi_downwash = atan2( in.AeroUVW(eV), v_induced - in.AeroUVW(eW));
+ calc_downwash_angles();
+ // ... and assign to inherited vFn and vMn members
+ // (for processing see FGForce::GetBodyForces).
vFn = body_forces(A_IC, B_IC);
- vMn = Transform() * body_moments(A_IC, B_IC);
+ vMn = Transform() * body_moments(A_IC, B_IC);
}
CalcRotorState();
if (! ExternalRPM) {
- Transmission->SetClutchCtrlNorm(ClutchCtrlNorm);
-
// the RPM values are handled inside Transmission
Transmission->Calculate(EnginePower, Torque, in.TotalDeltaT);
property_name = base_property_name + "/engine-rpm";
PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetEngineRPM );
- property_name = base_property_name + "/rotor-thrust-lbs"; // might be redundant - check!
- PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetThrust );
-
property_name = base_property_name + "/a0-rad";
PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetA0 );
property_name = base_property_name + "/phi-downwash-rad";
PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetPhiDW );
+ property_name = base_property_name + "/groundeffect-scale-norm";
+ PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetGroundEffectScaleNorm,
+ &FGRotor::SetGroundEffectScaleNorm );
+
switch (ControlMap) {
case eTailCtrl:
property_name = base_property_name + "/antitorque-ctrl-rad";
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-string FGRotor::GetThrusterLabels(int id, string delimeter)
+string FGRotor::GetThrusterLabels(int id, const string& delimeter)
{
ostringstream buf;
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-string FGRotor::GetThrusterValues(int id, string delimeter)
+string FGRotor::GetThrusterValues(int id, const string& delimeter)
{
ostringstream buf;
}
-} // namespace JSBSim
+} // namespace JSBSim
01/01/10 T.Kreitler test implementation
01/10/11 T.Kreitler changed to single rotor model
03/06/11 T.Kreitler added brake, clutch, and experimental free-wheeling-unit
+02/05/12 T.Kreitler brake, clutch, and FWU now in FGTransmission class
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
SENTRY
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
#include "FGThruster.h"
+#include "FGTransmission.h"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
DEFINITIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#define ID_ROTOR "$Id: FGRotor.h,v 1.12 2011/10/15 21:30:28 jentron Exp $"
+#define ID_ROTOR "$Id: FGRotor.h,v 1.14 2012/03/18 15:48:36 jentron Exp $"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FORWARD DECLARATIONS
<groundeffectexp> {number} </groundeffectexp>
<groundeffectshift unit="{LENGTH}"> {number} </groundeffectshift>
- <freewheelthresh> {number} </freewheelthresh>
</rotor>
// LENGTH means any of the supported units, same for ANGLE and MOMENT.
Experimental properties
\<groundeffectexp> - Exponent for ground effect approximation. Values usually range from 0.04
- for large rotors to 0.1 for smaller ones. As a rule of thumb the effect
+ for large rotors to 0.1 for smaller ones. As a rule of thumb the effect
vanishes at a height 2-3 times the rotor diameter.
formula used: exp ( - groundeffectexp * (height+groundeffectshift) )
Omitting or setting to 0.0 disables the effect calculation.
- \<groundeffectshift> - Further adjustment of ground effect, approx. hub height or slightly above.
+ \<groundeffectshift> - Further adjustment of ground effect, approx. hub height or slightly above
+ (This lessens the influence of the ground effect).
</pre>
-<h3>Notes:</h3>
+<h3>Notes:</h3>
<h4>- Controls -</h4>
The behavior of the rotor is controlled/influenced by following inputs.<ul>
<li> The power provided by the engine. This is handled by the regular engine controls.</li>
- <li> The collective control input. This is read from the <tt>fdm</tt> property
+ <li> The collective control input. This is read from the <tt>fdm</tt> property
<tt>propulsion/engine[x]/collective-ctrl-rad</tt>. See below for tail rotor</li>
<li> The lateral cyclic input. Read from
<tt>propulsion/engine[x]/lateral-ctrl-rad</tt>.</li>
<li> The longitudinal cyclic input. Read from
<tt>propulsion/engine[x]/longitudinal-ctrl-rad</tt>.</li>
- <li> The tail collective (aka antitorque, aka pedal) control input. Read from
- <tt>propulsion/engine[x]/antitorque-ctrl-rad</tt> or
+ <li> The tail rotor collective (aka antitorque, aka pedal) control input. Read from
+ <tt>propulsion/engine[x]/antitorque-ctrl-rad</tt> or
<tt>propulsion/engine[x]/tail-collective-ctrl-rad</tt>.</li>
</ul>
is linked to to the main (=first, =0) rotor, and specifing
<tt>\<controlmap\> TAIL \</controlmap\></tt> tells this rotor to read the
collective input from <tt>propulsion/engine[1]/antitorque-ctrl-rad</tt>
- (The TAIL-map ignores lateral and longitudinal input). The rotor needs to be
+ (The TAIL-map ignores lateral and longitudinal input). The rotor needs to be
attached to a dummy engine, e.g. an 1HP electrical engine.
A tandem rotor is setup analogous.
</ul>
+ <h4>- Scaling the ground effect -</h4>
+
+ The property <tt>propulsion/engine[x]/groundeffect-scale-norm</tt> allows fdm based
+ scaling of the ground effect influence. For instance the effect vanishes at speeds
+ above approx. 50kts, or one likes to land on a 'perforated' helipad.
+
<h4>- Development hints -</h4>
Setting <tt>\<ExternalRPM> -1 \</ExternalRPM></tt> the rotor's RPM is controlled by
when developing a FDM.
-<h3>References:</h3>
+<h3>References:</h3>
<dl>
<dt>/SH79/</dt><dd>Shaugnessy, J. D., Deaux, Thomas N., and Yenni, Kenneth R.,
- "Development and Validation of a Piloted Simulation of a
+ "Development and Validation of a Piloted Simulation of a
Helicopter and External Sling Load", NASA TP-1285, 1979.</dd>
<dt>/BA41/</dt><dd>Bailey,F.J.,Jr., "A Simplified Theoretical Method of Determining
the Characteristics of a Lifting Rotor in Forward Flight", NACA Rep.716, 1941.</dd>
<dt>/AM50/</dt><dd>Amer, Kenneth B.,"Theory of Helicopter Damping in Pitch or Roll and a
Comparison With Flight Measurements", NACA TN-2136, 1950.</dd>
<dt>/TA77/</dt><dd>Talbot, Peter D., Corliss, Lloyd D., "A Mathematical Force and Moment
- Model of a UH-1H Helicopter for Flight Dynamics Simulations", NASA TM-73,254, 1977.</dd>
+ Model of a UH-1H Helicopter for Flight Dynamics Simulations", NASA TM-73,254, 1977.</dd>
<dt>/GE49/</dt><dd>Gessow, Alfred, Amer, Kenneth B. "An Introduction to the Physical
- Aspects of Helicopter Stability", NACA TN-1982, 1949.</dd>
+ Aspects of Helicopter Stability", NACA TN-1982, 1949.</dd>
</dl>
@author Thomas Kreitler
- @version $Id: FGRotor.h,v 1.12 2011/10/15 21:30:28 jentron Exp $
+ @version $Id: FGRotor.h,v 1.14 2012/03/18 15:48:36 jentron Exp $
*/
CLASS DECLARATION
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-class FGTransmission : public FGJSBBase {
-
-public:
- FGTransmission(FGFDMExec *exec, int num);
- ~FGTransmission();
-
- void Calculate(double EnginePower, double ThrusterTorque, double dt);
-
- void SetMaxBrakePower(double x) {MaxBrakePower=x;}
- double GetMaxBrakePower() const {return MaxBrakePower;}
- void SetEngineFriction(double x) {EngineFriction=x;}
- double GetEngineFriction() const {return EngineFriction;}
- void SetEngineMoment(double x) {EngineMoment=x;}
- double GetEngineMoment() const {return EngineMoment;}
- void SetThrusterMoment(double x) {ThrusterMoment=x;}
- double GetThrusterMoment() const {return ThrusterMoment;}
-
- double GetFreeWheelTransmission() const {return FreeWheelTransmission;}
- double GetEngineRPM() {return EngineRPM;}
- double GetThrusterRPM() {return ThrusterRPM;}
-
- double GetBrakeCtrl() const {return BrakeCtrlNorm;}
- void SetBrakeCtrl(double x) {BrakeCtrlNorm=x;}
- void SetClutchCtrlNorm(double x) {ClutchCtrlNorm=x;}
-
-private:
- bool BindModel(int num);
- // void Debug(int from);
-
- inline double omega_to_rpm(double w) {return w * 9.54929658551372014613302580235;} // omega/(2.0*PI) * 60.0
- inline double rpm_to_omega(double r) {return r * .104719755119659774615421446109;} // (rpm/60.0)*2.0*PI
-
- Filter FreeWheelLag;
- double FreeWheelTransmission; // state, 0: free, 1:locked
-
- double ThrusterMoment;
- double EngineMoment; // estimated MOI of gear and engine, influences acceleration
- double EngineFriction; // estimated friction in gear and possibly engine
-
- double ClutchCtrlNorm; // also in FGThruster.h
- double BrakeCtrlNorm;
- double MaxBrakePower;
-
- double EngineRPM;
- double ThrusterRPM;
- FGPropertyManager* PropertyManager;
-
-};
-
-//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-
class FGRotor : public FGThruster {
enum eCtrlMapping {eMainCtrl=0, eTailCtrl, eTandemCtrl};
/// Destructor for FGRotor
~FGRotor();
- /** Returns the power required by the rotor. */
+ /// Returns the power required by the rotor.
double GetPowerRequired(void)const { return PowerRequired; }
- /** Returns the scalar thrust of the rotor, and adjusts the RPM value. */
+ /// Returns the scalar thrust of the rotor, and adjusts the RPM value.
double Calculate(double EnginePower);
/// Retrieves the torque
double GetTorque(void) const { return Torque; }
- /// Downwash angle - currently only valid for a rotor that spins horizontally
+ /// Downwash angle - positive values point forward (given a horizontal spinning rotor)
double GetThetaDW(void) const { return theta_downwash; }
- /// Downwash angle - currently only valid for a rotor that spins horizontally
+ /// Downwash angle - positive values point leftward (given a horizontal spinning rotor)
double GetPhiDW(void) const { return phi_downwash; }
+ /// Retrieves the ground effect scaling factor.
+ double GetGroundEffectScaleNorm(void) const { return GroundEffectScaleNorm; }
+ /// Sets the ground effect scaling factor.
+ void SetGroundEffectScaleNorm(double g) { GroundEffectScaleNorm = g; }
+
/// Retrieves the collective control input in radians.
double GetCollectiveCtrl(void) const { return CollectiveCtrl; }
/// Retrieves the lateral control input in radians.
void SetLongitudinalCtrl(double c) { LongitudinalCtrl = c; }
// Stubs. Only main rotor RPM is returned
- string GetThrusterLabels(int id, string delimeter);
- string GetThrusterValues(int id, string delimeter);
+ string GetThrusterLabels(int id, const string& delimeter);
+ string GetThrusterValues(int id, const string& delimeter);
private:
double ConfigValue( Element* e, const string& ename, double default_val=0.0,
bool tell=false);
- void Configure(Element* rotor_element);
+ double Configure(Element* rotor_element);
void CalcRotorState(void);
void calc_flapping_angles(double theta_0, const FGColumnVector3 &pqr_fus_w);
void calc_drag_and_side_forces(double theta_0);
void calc_torque(double theta_0);
+ void calc_downwash_angles();
// transformations
FGColumnVector3 hub_vel_body2ca( const FGColumnVector3 &uvw, const FGColumnVector3 &pqr,
FGPropertyManager* ExtRPMsource;
double SourceGearRatio;
+ // 'real' rotor parameters
double BladeChord;
double LiftCurveSlope;
double BladeTwist;
double InflowLag;
double TipLossB;
+ // groundeffect
double GroundEffectExp;
double GroundEffectShift;
+ double GroundEffectScaleNorm;
// derived parameters
double LockNumberByRho;
double lambda; // inflow ratio
double mu; // tip-speed ratio
double nu; // induced inflow ratio
- double v_induced; // induced velocity, always positive [ft/s]
+ double v_induced; // induced velocity, usually positive [ft/s]
double theta_downwash;
double phi_downwash;
namespace JSBSim {
-static const char *IdSrc = "$Id: FGThruster.cpp,v 1.14 2011/03/10 01:35:25 dpculp Exp $";
+static const char *IdSrc = "$Id: FGThruster.cpp,v 1.16 2012/03/18 15:48:36 jentron Exp $";
static const char *IdHdr = ID_THRUSTER;
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
GearRatio = 1.0;
ReverserAngle = 0.0;
- ClutchCtrlNorm = 1.0;
EngineNum = num;
PropertyManager = FDMExec->GetPropertyManager();
&FGThruster::SetReverserAngle);
}
- if (el->GetName() == "rotor") // At this time only a rotor can have a clutch.
- {
- property_name = base_property_name + "/clutch-ctrl-norm";
- PropertyManager->Tie( property_name.c_str(), (FGThruster *)this, &FGThruster::GetClutchCtrl,
- &FGThruster::SetClutchCtrl);
- }
-
Debug(0);
}
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-string FGThruster::GetThrusterLabels(int id, string delimeter)
+string FGThruster::GetThrusterLabels(int id, const string& delimeter)
{
std::ostringstream buf;
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-string FGThruster::GetThrusterValues(int id, string delimeter)
+string FGThruster::GetThrusterValues(int id, const string& delimeter)
{
std::ostringstream buf;
DEFINITIONS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
-#define ID_THRUSTER "$Id: FGThruster.h,v 1.18 2011/09/24 14:26:46 jentron Exp $"
+#define ID_THRUSTER "$Id: FGThruster.h,v 1.20 2012/03/18 15:48:36 jentron Exp $"
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
FORWARD DECLARATIONS
1.57 (pi/2) results in no thrust at all.
@author Jon Berndt
- @version $Id: FGThruster.h,v 1.18 2011/09/24 14:26:46 jentron Exp $
+ @version $Id: FGThruster.h,v 1.20 2012/03/18 15:48:36 jentron Exp $
*/
/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
string GetName(void) {return Name;}
void SetReverserAngle(double angle) {ReverserAngle = angle;}
double GetReverserAngle(void) const {return ReverserAngle;}
- double GetClutchCtrl(void) const { return ClutchCtrlNorm; }
- void SetClutchCtrl(double c) { ClutchCtrlNorm = c; }
virtual double GetRPM(void) const { return 0.0; };
virtual double GetEngineRPM(void) const { return 0.0; };
double GetGearRatio(void) {return GearRatio; }
- virtual string GetThrusterLabels(int id, string delimeter);
- virtual string GetThrusterValues(int id, string delimeter);
+ virtual string GetThrusterLabels(int id, const string& delimeter);
+ virtual string GetThrusterValues(int id, const string& delimeter);
struct Inputs {
double TotalDeltaT;
double GearRatio;
double ThrustCoeff;
double ReverserAngle;
- double ClutchCtrlNorm;
int EngineNum;
FGPropertyManager* PropertyManager;
virtual void Debug(int from);
--- /dev/null
+/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+ JSBSim
+ Author: Jon S. Berndt
+ Date started: 08/24/00
+ ------------- Copyright (C) 2000 Jon S. Berndt (jon@jsbsim.org) -------------
+
+ Module: FGTransmission.cpp
+ Author: T.Kreitler
+ Date started: 02/05/12
+
+ ------------- Copyright (C) 2012 T. Kreitler (t.kreitler@web.de) -------------
+
+ This program is free software; you can redistribute it and/or modify it under
+ the terms of the GNU Lesser General Public License as published by the Free Software
+ Foundation; either version 2 of the License, or (at your option) any later
+ version.
+
+ This program is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+ FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
+ details.
+
+ You should have received a copy of the GNU Lesser General Public License along with
+ this program; if not, write to the Free Software Foundation, Inc., 59 Temple
+ Place - Suite 330, Boston, MA 02111-1307, USA.
+
+ Further information about the GNU Lesser General Public License can also be found on
+ the world wide web at http://www.gnu.org.
+
+FUNCTIONAL DESCRIPTION
+--------------------------------------------------------------------------------
+
+HISTORY
+--------------------------------------------------------------------------------
+02/05/12 T.Kreitler Created
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+INCLUDES
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
+
+
+#include "FGTransmission.h"
+
+using std::cout;
+using std::endl;
+
+namespace JSBSim {
+
+static const char *IdSrc = "$Id: FGTransmission.cpp,v 1.1 2012/02/25 14:37:02 jentron Exp $";
+static const char *IdHdr = ID_TRANSMISSION;
+
+/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+CLASS IMPLEMENTATION
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
+
+FGTransmission::FGTransmission(FGFDMExec *exec, int num, double dt) :
+ FreeWheelTransmission(1.0),
+ ThrusterMoment(1.0), EngineMoment(1.0), EngineFriction(0.0),
+ ClutchCtrlNorm(1.0), BrakeCtrlNorm(0.0), MaxBrakePower(0.0),
+ EngineRPM(0.0), ThrusterRPM(0.0)
+{
+ PropertyManager = exec->GetPropertyManager();
+ FreeWheelLag = Filter(200.0,dt); // avoid too abrupt changes in transmission
+ BindModel(num);
+}
+
+//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+FGTransmission::~FGTransmission(){
+}
+
+//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+// basically P = Q*w and Q_Engine + (-Q_Rotor) = J * dw/dt, J = Moment
+//
+void FGTransmission::Calculate(double EnginePower, double ThrusterTorque, double dt) {
+
+ double coupling = 1.0, coupling_sq = 1.0;
+ double fw_mult = 1.0;
+
+ double d_omega = 0.0, engine_d_omega = 0.0, thruster_d_omega = 0.0; // relative changes
+
+ double engine_omega = rpm_to_omega(EngineRPM);
+ double safe_engine_omega = engine_omega < 1e-1 ? 1e-1 : engine_omega;
+ double engine_torque = EnginePower / safe_engine_omega;
+
+ double thruster_omega = rpm_to_omega(ThrusterRPM);
+ double safe_thruster_omega = thruster_omega < 1e-1 ? 1e-1 : thruster_omega;
+
+ engine_torque -= EngineFriction / safe_engine_omega;
+ ThrusterTorque += Constrain(0.0, BrakeCtrlNorm, 1.0) * MaxBrakePower / safe_thruster_omega;
+
+ // would the FWU release ?
+ engine_d_omega = engine_torque/EngineMoment * dt;
+ thruster_d_omega = - ThrusterTorque/ThrusterMoment * dt;
+
+ if ( thruster_omega+thruster_d_omega > engine_omega+engine_d_omega ) {
+ // don't drive the engine
+ FreeWheelTransmission = 0.0;
+ } else {
+ FreeWheelTransmission = 1.0;
+ }
+
+ fw_mult = FreeWheelLag.execute(FreeWheelTransmission);
+ coupling = fw_mult * Constrain(0.0, ClutchCtrlNorm, 1.0);
+
+ if (coupling < 0.999999) { // are the separate calculations needed ?
+ // assume linear transfer
+ engine_d_omega =
+ (engine_torque - ThrusterTorque*coupling)/(ThrusterMoment*coupling + EngineMoment) * dt;
+ thruster_d_omega =
+ (engine_torque*coupling - ThrusterTorque)/(ThrusterMoment + EngineMoment*coupling) * dt;
+
+ EngineRPM += omega_to_rpm(engine_d_omega);
+ ThrusterRPM += omega_to_rpm(thruster_d_omega);
+
+ // simulate transit to static friction
+ coupling_sq = coupling*coupling;
+ EngineRPM = (1.0-coupling_sq) * EngineRPM + coupling_sq * 0.02 * (49.0*EngineRPM + ThrusterRPM);
+ ThrusterRPM = (1.0-coupling_sq) * ThrusterRPM + coupling_sq * 0.02 * (EngineRPM + 49.0*ThrusterRPM);
+
+ // enforce equal rpm
+ if ( fabs(EngineRPM-ThrusterRPM) < 1e-3 ) {
+ EngineRPM = ThrusterRPM = 0.5 * (EngineRPM+ThrusterRPM);
+ }
+ } else {
+ d_omega = (engine_torque - ThrusterTorque)/(ThrusterMoment + EngineMoment) * dt;
+ EngineRPM = ThrusterRPM += omega_to_rpm(d_omega);
+ }
+
+ // nothing will turn backward
+ if (EngineRPM < 0.0 ) EngineRPM = 0.0;
+ if (ThrusterRPM < 0.0 ) ThrusterRPM = 0.0;
+
+}
+
+//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+bool FGTransmission::BindModel(int num)
+{
+ string property_name, base_property_name;
+ base_property_name = CreateIndexedPropertyName("propulsion/engine", num);
+
+ property_name = base_property_name + "/brake-ctrl-norm";
+ PropertyManager->Tie( property_name.c_str(), this,
+ &FGTransmission::GetBrakeCtrlNorm, &FGTransmission::SetBrakeCtrlNorm);
+
+ property_name = base_property_name + "/clutch-ctrl-norm";
+ PropertyManager->Tie( property_name.c_str(), this,
+ &FGTransmission::GetClutchCtrlNorm, &FGTransmission::SetClutchCtrlNorm);
+
+ property_name = base_property_name + "/free-wheel-transmission";
+ PropertyManager->Tie( property_name.c_str(), this,
+ &FGTransmission::GetFreeWheelTransmission);
+
+ return true;
+}
+
+//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+// The bitmasked value choices are as follows:
+// unset: In this case (the default) JSBSim would only print
+// out the normally expected messages, essentially echoing
+// the config files as they are read. If the environment
+// variable is not set, debug_lvl is set to 1 internally
+// 0: This requests JSBSim not to output any messages
+// whatsoever.
+// 1: This value explicity requests the normal JSBSim
+// startup messages
+// 2: This value asks for a message to be printed out when
+// a class is instantiated
+// 4: When this value is set, a message is displayed when a
+// FGModel object executes its Run() method
+// 8: When this value is set, various runtime state variables
+// are printed out periodically
+// 16: When set various parameters are sanity checked and
+// a message is printed out when they go out of bounds
+
+void FGTransmission::Debug(int from)
+{
+ if (debug_lvl <= 0) return;
+
+ if (debug_lvl & 1) { // Standard console startup message output
+ if (from == 0) { // Constructor
+ }
+ }
+ if (debug_lvl & 2 ) { // Instantiation/Destruction notification
+ if (from == 0) cout << "Instantiated: FGTransmission" << endl;
+ if (from == 1) cout << "Destroyed: FGTransmission" << endl;
+ }
+ if (debug_lvl & 4 ) { // Run() method entry print for FGModel-derived objects
+ }
+ if (debug_lvl & 8 ) { // Runtime state variables
+ }
+ if (debug_lvl & 16) { // Sanity checking
+ }
+ if (debug_lvl & 64) {
+ if (from == 0) { // Constructor
+ cout << IdSrc << endl;
+ cout << IdHdr << endl;
+ }
+ }
+}
+
+}
+
--- /dev/null
+/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+ JSBSim
+ Author: Jon S. Berndt
+ Date started: 08/24/00
+ ------------- Copyright (C) 2000 Jon S. Berndt (jon@jsbsim.org) -------------
+
+ Header: FGTransmission.h
+ Author: T.Kreitler
+ Date started: 02/05/12
+
+ ------------- Copyright (C) 2012 T. Kreitler (t.kreitler@web.de) -------------
+
+ This program is free software; you can redistribute it and/or modify it under
+ the terms of the GNU Lesser General Public License as published by the Free Software
+ Foundation; either version 2 of the License, or (at your option) any later
+ version.
+
+ This program is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+ FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
+ details.
+
+ You should have received a copy of the GNU Lesser General Public License along with
+ this program; if not, write to the Free Software Foundation, Inc., 59 Temple
+ Place - Suite 330, Boston, MA 02111-1307, USA.
+
+ Further information about the GNU Lesser General Public License can also be found on
+ the world wide web at http://www.gnu.org.
+
+HISTORY
+--------------------------------------------------------------------------------
+02/05/12 T.Kreitler Created
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+SENTRY
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
+
+#ifndef FGTRANSMISSION_H
+#define FGTRANSMISSION_H
+
+/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+INCLUDES
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
+
+#include "FGJSBBase.h"
+#include "FGFDMExec.h"
+#include "input_output/FGPropertyManager.h"
+
+/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+DEFINITIONS
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
+
+#define ID_TRANSMISSION "$Id: FGTransmission.h,v 1.1 2012/02/25 14:37:02 jentron Exp $"
+
+/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+FORWARD DECLARATIONS
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
+
+namespace JSBSim {
+
+/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+CLASS DOCUMENTATION
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
+
+/** Utility class that handles power transmission in conjunction with FGRotor.
+
+ This class provides brake, clutch and free-wheel-unit (FWU) functionality
+ for the rotor model. Also it is responsible for the RPM calculations.
+
+ When the engine is off the brake could be used to slow/hold down a spinning
+ rotor. The maximum brake power is defined in the rotors' config file.
+ (Right now there is no checking if the input is in the [0..1] range.)
+
+ The clutch operation is based on a heuristic approach. In the intermediate
+ state the transfer is proportional to the clutch position. But equal RPM
+ values are enforced on the thruster and rotor sides when approaching the
+ closed state.
+
+ The FWU inhibits that the rotor is driving the engine. To do so, the code
+ just predicts the upcoming FWU state based on current torque conditions.
+
+ Some engines won't work properly when the clutch is open. To keep them
+ controllable some load must be provided on the engine side (EngineFriction,
+ aka gear-loss). See the notes under 'Engine issues' in FGRotor.
+
+<h3>Property tree</h3>
+
+ The following properties are created (with x = your thruster number):
+ <pre>
+ propulsion/engine[x]/brake-ctrl-norm
+ propulsion/engine[x]/free-wheel-transmission
+ propulsion/engine[x]/clutch-ctrl-norm
+ </pre>
+
+<h3>Notes</h3>
+
+ <ul>
+ <li> EngineFriction is assumed constant, so better orientate at low RPM
+ values, because piston and turboprop engines don't 'like' high
+ load at startup.</li>
+ <li> The model doesn't support backward operation.</li>
+ <li> And even worse, the torque calculations silently assume a minimal
+ RPM value of approx. 1.</li>
+ </ul>
+
+*/
+
+
+/*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+CLASS DECLARATION
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
+
+class FGTransmission : public FGJSBBase {
+
+public:
+ /** Constructor for FGTransmission.
+ @param exec a pointer to the main executive object
+ @param num the number of the thruster that uses this object
+ @param dt simulation delta T */
+ FGTransmission(FGFDMExec *exec, int num, double dt);
+
+ /// Destructor for FGTransmission
+ ~FGTransmission();
+
+ void Calculate(double EnginePower, double ThrusterTorque, double dt);
+
+ void SetMaxBrakePower(double x) {MaxBrakePower=x;}
+ double GetMaxBrakePower() const {return MaxBrakePower;}
+ void SetEngineFriction(double x) {EngineFriction=x;}
+ double GetEngineFriction() const {return EngineFriction;}
+ void SetEngineMoment(double x) {EngineMoment=x;}
+ double GetEngineMoment() const {return EngineMoment;}
+ void SetThrusterMoment(double x) {ThrusterMoment=x;}
+ double GetThrusterMoment() const {return ThrusterMoment;}
+
+ double GetFreeWheelTransmission() const {return FreeWheelTransmission;}
+ void SetEngineRPM(double x) {EngineRPM=x;}
+ double GetEngineRPM() {return EngineRPM;}
+ void SetThrusterRPM(double x) {ThrusterRPM=x;}
+ double GetThrusterRPM() {return ThrusterRPM;}
+
+ double GetBrakeCtrlNorm() const {return BrakeCtrlNorm;}
+ void SetBrakeCtrlNorm(double x) {BrakeCtrlNorm=x;}
+ double GetClutchCtrlNorm() const {return ClutchCtrlNorm;}
+ void SetClutchCtrlNorm(double x) {ClutchCtrlNorm=x;}
+
+private:
+ bool BindModel(int num);
+ void Debug(int from);
+
+ inline double omega_to_rpm(double w) {
+ return w * 9.54929658551372014613302580235; // omega/(2.0*PI) * 60.0
+ }
+ inline double rpm_to_omega(double r) {
+ return r * 0.104719755119659774615421446109; // (rpm/60.0)*2.0*PI
+ }
+
+ Filter FreeWheelLag;
+ double FreeWheelTransmission; // state, 0: free, 1:locked
+
+ double ThrusterMoment;
+ double EngineMoment; // estimated MOI of gear and engine, influences acceleration
+ double EngineFriction; // estimated friction in gear and possibly engine
+
+ double ClutchCtrlNorm;
+ double BrakeCtrlNorm;
+ double MaxBrakePower;
+
+ double EngineRPM;
+ double ThrusterRPM;
+ FGPropertyManager* PropertyManager;
+
+};
+
+}
+//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+#endif
+
projects / flightgear.git / commitdiff