diff --git a/src/FDM/JSBSim/FGFDMExec.cpp b/src/FDM/JSBSim/FGFDMExec.cpp index 2d5d409fa..25f849a9d 100644 --- a/src/FDM/JSBSim/FGFDMExec.cpp +++ b/src/FDM/JSBSim/FGFDMExec.cpp @@ -572,6 +572,9 @@ void FGFDMExec::BuildPropertyCatalog(struct PropertyCatalogStructure* pcs) sprintf(int_buf, "[%d]", node_idx); if (node_idx != 0) pcsNew->base_string += string(int_buf); if (pcs->node->getChild(i)->nChildren() == 0) { + if (pcsNew->base_string.substr(0,11) == string("/fdm/jsbsim")) { + pcsNew->base_string = pcsNew->base_string.erase(0,12); + } PropertyCatalog.push_back(pcsNew->base_string); } else { pcsNew->node = (FGPropertyManager*)pcs->node->getChild(i); diff --git a/src/FDM/JSBSim/JSBSim.cxx b/src/FDM/JSBSim/JSBSim.cxx index f63e106ad..a13821d82 100644 --- a/src/FDM/JSBSim/JSBSim.cxx +++ b/src/FDM/JSBSim/JSBSim.cxx @@ -282,6 +282,10 @@ FGJSBsim::FGJSBsim( double dt ) fgGetDouble("/fdm/jsbsim/systems/hook/tailhook-offset-y-in", 0), fgGetDouble("/fdm/jsbsim/systems/hook/tailhook-offset-z-in", -16)); + // Untie the write-state-file property to avoid creating an initfile.xml + // file on each reset. + fgGetNode("/fdm/jsbsim/simulation/write-state-file")->untie(); + crashed = false; } diff --git a/src/FDM/JSBSim/initialization/FGInitialCondition.cpp b/src/FDM/JSBSim/initialization/FGInitialCondition.cpp index eedec9345..2fbbae224 100644 --- a/src/FDM/JSBSim/initialization/FGInitialCondition.cpp +++ b/src/FDM/JSBSim/initialization/FGInitialCondition.cpp @@ -156,7 +156,13 @@ void FGInitialCondition::InitializeIC(void) void FGInitialCondition::WriteStateFile(int num) { - string filename = fdmex->GetFullAircraftPath() + "/" + "initfile.xml"; + string filename = fdmex->GetFullAircraftPath(); + + if (filename.empty()) + filename = "initfile.xml"; + else + filename.append("/initfile.xml"); + ofstream outfile(filename.c_str()); FGPropagate* Propagate = fdmex->GetPropagate(); @@ -173,11 +179,10 @@ void FGInitialCondition::WriteStateFile(int num) outfile << " " << Propagate->GetLatitudeDeg() << " " << endl; outfile << " " << Propagate->Geth() << " " << endl; outfile << "" << endl; + outfile.close(); } else { - cerr << "Could not open and/or write the state to the initial conditions file." << endl; + cerr << "Could not open and/or write the state to the initial conditions file: " << filename << endl; } - - outfile.close(); } //****************************************************************************** diff --git a/src/FDM/JSBSim/input_output/FGScript.cpp b/src/FDM/JSBSim/input_output/FGScript.cpp index e0cddc34b..9ea8a9186 100755 --- a/src/FDM/JSBSim/input_output/FGScript.cpp +++ b/src/FDM/JSBSim/input_output/FGScript.cpp @@ -46,7 +46,6 @@ INCLUDES #include #include -#include namespace JSBSim { @@ -204,12 +203,17 @@ bool FGScript::LoadScript( string script ) // Retrieve the event name if given newEvent->Name = event_element->GetAttributeValue("name"); - // Is this event persistent? That is, does it execute repeatedly as long as the - // condition is true, or does it execute as a one-shot event, only? + // Is this event persistent? That is, does it execute every time the + // condition triggers to true, or does it execute as a one-shot event, only? if (event_element->GetAttributeValue("persistent") == string("true")) { newEvent->Persistent = true; } + // Does this event execute continuously when triggered to true? + if (event_element->GetAttributeValue("continuous") == string("true")) { + newEvent->Continuous = true; + } + // Process the conditions condition_element = event_element->FindElement("condition"); if (condition_element != 0) { @@ -304,7 +308,6 @@ bool FGScript::LoadScript( string script ) bool FGScript::RunScript(void) { - vector ::iterator iEvent = Events.begin(); unsigned i, j; unsigned event_ctr = 0; @@ -314,90 +317,98 @@ bool FGScript::RunScript(void) if (currentTime > EndTime) return false; //Script done! // Iterate over all events. - while (iEvent < Events.end()) { - iEvent->PrevTriggered = iEvent->Triggered; + for (unsigned int ev_ctr=0; ev_ctr < Events.size(); ev_ctr++) { // Determine whether the set of conditional tests for this condition equate - // to true and should cause the event to execute. - if (iEvent->Condition->Evaluate()) { - if (!iEvent->Triggered) { + // to true and should cause the event to execute. If the conditions evaluate + // to true, then the event is triggered. If the event is not persistent, + // then this trigger will remain set true. If the event is persistent, + // the trigger will reset to false when the condition evaluates to false. + if (Events[ev_ctr].Condition->Evaluate()) { + if (!Events[ev_ctr].Triggered) { // The conditions are true, do the setting of the desired Event parameters - for (i=0; iSetValue.size(); i++) { - iEvent->OriginalValue[i] = iEvent->SetParam[i]->getDoubleValue(); - if (iEvent->Functions[i] != 0) { // Parameter should be set to a function value - iEvent->SetValue[i] = iEvent->Functions[i]->GetValue(); + for (i=0; igetDoubleValue(); + if (Events[ev_ctr].Functions[i] != 0) { // Parameter should be set to a function value + Events[ev_ctr].SetValue[i] = Events[ev_ctr].Functions[i]->GetValue(); } - switch (iEvent->Type[i]) { + switch (Events[ev_ctr].Type[i]) { case FG_VALUE: case FG_BOOL: - iEvent->newValue[i] = iEvent->SetValue[i]; + Events[ev_ctr].newValue[i] = Events[ev_ctr].SetValue[i]; break; case FG_DELTA: - iEvent->newValue[i] = iEvent->OriginalValue[i] + iEvent->SetValue[i]; + Events[ev_ctr].newValue[i] = Events[ev_ctr].OriginalValue[i] + Events[ev_ctr].SetValue[i]; break; default: cerr << "Invalid Type specified" << endl; break; } - iEvent->StartTime = currentTime + iEvent->Delay; - iEvent->ValueSpan[i] = iEvent->newValue[i] - iEvent->OriginalValue[i]; - iEvent->Transiting[i] = true; + Events[ev_ctr].StartTime = currentTime + Events[ev_ctr].Delay; + Events[ev_ctr].ValueSpan[i] = Events[ev_ctr].newValue[i] - Events[ev_ctr].OriginalValue[i]; + Events[ev_ctr].Transiting[i] = true; } } - iEvent->Triggered = true; - } else if (iEvent->Persistent) { - iEvent->Triggered = false; // Reset the trigger for persistent events - iEvent->Notified = false; // Also reset the notification flag + Events[ev_ctr].Triggered = true; + + } else if (Events[ev_ctr].Persistent) { // If the event is persistent, reset the trigger. + + Events[ev_ctr].Triggered = false; // Reset the trigger for persistent events + Events[ev_ctr].Notified = false; // Also reset the notification flag } - if ((currentTime >= iEvent->StartTime) && iEvent->Triggered) { + if ((currentTime >= Events[ev_ctr].StartTime) && Events[ev_ctr].Triggered) { - for (i=0; iSetValue.size(); i++) { - if (iEvent->Transiting[i]) { - iEvent->TimeSpan = currentTime - iEvent->StartTime; - if (iEvent->Functions[i] == 0) { - switch (iEvent->Action[i]) { - case FG_RAMP: - if (iEvent->TimeSpan <= iEvent->TC[i]) { - newSetValue = iEvent->TimeSpan/iEvent->TC[i] * iEvent->ValueSpan[i] + iEvent->OriginalValue[i]; - } else { - newSetValue = iEvent->newValue[i]; - iEvent->Transiting[i] = false; - } - break; - case FG_STEP: - newSetValue = iEvent->newValue[i]; - iEvent->Transiting[i] = false; - break; - case FG_EXP: - newSetValue = (1 - exp( -iEvent->TimeSpan/iEvent->TC[i] )) * iEvent->ValueSpan[i] + iEvent->OriginalValue[i]; - break; - default: - cerr << "Invalid Action specified" << endl; - break; + for (i=0; iFunctions[i]->GetValue(); + break; + case FG_STEP: + newSetValue = Events[ev_ctr].newValue[i]; + + // If this is not a continuous event, reset the transiting flag. + // Otherwise, it is known that the event is a continuous event. + // Furthermore, if the event is to be determined by a function, + // then the function will be continuously calculated. + if (Events[ev_ctr].Continuous != true) + Events[ev_ctr].Transiting[i] = false; + else if (Events[ev_ctr].Functions[i] != 0) + newSetValue = Events[ev_ctr].Functions[i]->GetValue(); + + break; + case FG_EXP: + newSetValue = (1 - exp( -Events[ev_ctr].TimeSpan/Events[ev_ctr].TC[i] )) * Events[ev_ctr].ValueSpan[i] + Events[ev_ctr].OriginalValue[i]; + break; + default: + cerr << "Invalid Action specified" << endl; + break; } - iEvent->SetParam[i]->setDoubleValue(newSetValue); + Events[ev_ctr].SetParam[i]->setDoubleValue(newSetValue); } } // Print notification values after setting them - if (iEvent->Notify && !iEvent->Notified) { - cout << endl << " Event " << event_ctr << " (" << iEvent->Name << ")" + if (Events[ev_ctr].Notify && !Events[ev_ctr].Notified) { + cout << endl << " Event " << event_ctr << " (" << Events[ev_ctr].Name << ")" << " executed at time: " << currentTime << endl; - for (j=0; jNotifyProperties.size();j++) { - cout << " " << iEvent->NotifyProperties[j]->GetName() - << " = " << iEvent->NotifyProperties[j]->getDoubleValue() << endl; + for (j=0; jGetName() + << " = " << Events[ev_ctr].NotifyProperties[j]->getDoubleValue() << endl; } cout << endl; - iEvent->Notified = true; + Events[ev_ctr].Notified = true; } } - iEvent++; event_ctr++; } return true; diff --git a/src/FDM/JSBSim/input_output/FGScript.h b/src/FDM/JSBSim/input_output/FGScript.h index f1a09bf1f..fdb0a97ff 100644 --- a/src/FDM/JSBSim/input_output/FGScript.h +++ b/src/FDM/JSBSim/input_output/FGScript.h @@ -70,13 +70,14 @@ CLASS DOCUMENTATION format. A test condition (or conditions) can be set up in an event in a script and when the condition evaluates to true, the specified action[s] is/are taken. An event can be persistent, - meaning that at all times when the test condition evaluates to true - the specified set actions take place. When the set of - tests evaluates to true for a given + meaning that at every time the test condition first evaluates to true + (toggling from false to true) then the specified set actions take + place. An event can also be defined to execute or evaluate continuously + while the condition is true. When the set of tests evaluates to true for a given condition, an item may be set to another value. This value may be a value, or a delta value, and the change from the - current value to the new value can be either via a step function, - a ramp, or an exponential approach. The speed of a ramp or + current value to the new value can be either via a step action, + a ramp, or an exponential approach. The speed of a ramp or exponential approach is specified via the time constant. Here is an example illustrating the format of the script file: @@ -204,8 +205,8 @@ private: struct event { FGCondition *Condition; bool Persistent; + bool Continuous; bool Triggered; - bool PrevTriggered; bool Notify; bool Notified; double Delay; @@ -226,8 +227,8 @@ private: event() { Triggered = false; - PrevTriggered = false; Persistent = false; + Continuous = false; Delay = 0.0; Notify = Notified = false; Name = ""; @@ -237,7 +238,6 @@ private: void reset(void) { Triggered = false; - PrevTriggered = false; Notified = false; StartTime = 0.0; } diff --git a/src/FDM/JSBSim/math/FGColumnVector3.cpp b/src/FDM/JSBSim/math/FGColumnVector3.cpp index 0bd5c011f..0bdeadaba 100644 --- a/src/FDM/JSBSim/math/FGColumnVector3.cpp +++ b/src/FDM/JSBSim/math/FGColumnVector3.cpp @@ -52,7 +52,7 @@ CLASS IMPLEMENTATION FGColumnVector3::FGColumnVector3(void) { data[0] = data[1] = data[2] = 0.0; - Debug(0); + // Debug(0); } //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% @@ -60,7 +60,7 @@ FGColumnVector3::FGColumnVector3(void) string FGColumnVector3::Dump(string delimeter) const { char buffer[256]; - sprintf(buffer, "%10.3e%s%10.3e%s%10.3e", Entry(1), delimeter.c_str(), Entry(2), delimeter.c_str(), Entry(3)); + sprintf(buffer, "%13.6e%s%13.6e%s%13.6e", Entry(1), delimeter.c_str(), Entry(2), delimeter.c_str(), Entry(3)); return string(buffer); } diff --git a/src/FDM/JSBSim/math/FGMatrix33.cpp b/src/FDM/JSBSim/math/FGMatrix33.cpp index 9b4ce1231..91976f100 100644 --- a/src/FDM/JSBSim/math/FGMatrix33.cpp +++ b/src/FDM/JSBSim/math/FGMatrix33.cpp @@ -56,7 +56,7 @@ FGMatrix33::FGMatrix33(void) data[0] = data[1] = data[2] = data[3] = data[4] = data[5] = data[6] = data[7] = data[8] = 0.0; - Debug(0); + // Debug(0); } //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% diff --git a/src/FDM/JSBSim/models/FGAircraft.cpp b/src/FDM/JSBSim/models/FGAircraft.cpp index 5dd8612c6..463e902a2 100644 --- a/src/FDM/JSBSim/models/FGAircraft.cpp +++ b/src/FDM/JSBSim/models/FGAircraft.cpp @@ -140,7 +140,7 @@ bool FGAircraft::Run(void) //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -double FGAircraft::GetNlf(void) +double FGAircraft::GetNlf(void) const { return -1*Aerodynamics->GetvFw(3)/MassBalance->GetWeight(); } @@ -216,13 +216,6 @@ void FGAircraft::bind(void) PropertyManager->Tie("metrics/lv-norm", this, &FGAircraft::Getlbarv); PropertyManager->Tie("metrics/vbarh-norm", this, &FGAircraft::Getvbarh); PropertyManager->Tie("metrics/vbarv-norm", this, &FGAircraft::Getvbarv); - PropertyManager->Tie("forces/hold-down", this, &FGAircraft::GetHoldDown, &FGAircraft::SetHoldDown); - PropertyManager->Tie("moments/l-total-lbsft", this, eL, (PMF)&FGAircraft::GetMoments); - PropertyManager->Tie("moments/m-total-lbsft", this, eM, (PMF)&FGAircraft::GetMoments); - PropertyManager->Tie("moments/n-total-lbsft", this, eN, (PMF)&FGAircraft::GetMoments); - PropertyManager->Tie("forces/fbx-total-lbs", this, eX, (PMF)&FGAircraft::GetForces); - PropertyManager->Tie("forces/fby-total-lbs", this, eY, (PMF)&FGAircraft::GetForces); - PropertyManager->Tie("forces/fbz-total-lbs", this, eZ, (PMF)&FGAircraft::GetForces); PropertyManager->Tie("metrics/aero-rp-x-in", this, eX, (PMF)&FGAircraft::GetXYZrp); PropertyManager->Tie("metrics/aero-rp-y-in", this, eY, (PMF)&FGAircraft::GetXYZrp); PropertyManager->Tie("metrics/aero-rp-z-in", this, eZ, (PMF)&FGAircraft::GetXYZrp); @@ -232,6 +225,14 @@ void FGAircraft::bind(void) PropertyManager->Tie("metrics/visualrefpoint-x-in", this, eX, (PMF)&FGAircraft::GetXYZvrp); PropertyManager->Tie("metrics/visualrefpoint-y-in", this, eY, (PMF)&FGAircraft::GetXYZvrp); PropertyManager->Tie("metrics/visualrefpoint-z-in", this, eZ, (PMF)&FGAircraft::GetXYZvrp); + PropertyManager->Tie("forces/fbx-total-lbs", this, eX, (PMF)&FGAircraft::GetForces); + PropertyManager->Tie("forces/fby-total-lbs", this, eY, (PMF)&FGAircraft::GetForces); + PropertyManager->Tie("forces/fbz-total-lbs", this, eZ, (PMF)&FGAircraft::GetForces); + PropertyManager->Tie("forces/load-factor", this, &FGAircraft::GetNlf); + PropertyManager->Tie("forces/hold-down", this, &FGAircraft::GetHoldDown, &FGAircraft::SetHoldDown); + PropertyManager->Tie("moments/l-total-lbsft", this, eL, (PMF)&FGAircraft::GetMoments); + PropertyManager->Tie("moments/m-total-lbsft", this, eM, (PMF)&FGAircraft::GetMoments); + PropertyManager->Tie("moments/n-total-lbsft", this, eN, (PMF)&FGAircraft::GetMoments); } //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% diff --git a/src/FDM/JSBSim/models/FGAircraft.h b/src/FDM/JSBSim/models/FGAircraft.h index f833c3a79..616da6cf0 100644 --- a/src/FDM/JSBSim/models/FGAircraft.h +++ b/src/FDM/JSBSim/models/FGAircraft.h @@ -171,7 +171,7 @@ public: void SetWingArea(double S) {WingArea = S;} - double GetNlf(void); + double GetNlf(void) const; inline FGColumnVector3& GetNwcg(void) { return vNwcg; } diff --git a/src/FDM/JSBSim/models/FGAtmosphere.cpp b/src/FDM/JSBSim/models/FGAtmosphere.cpp index 88ffbd5e0..215b3b8a1 100644 --- a/src/FDM/JSBSim/models/FGAtmosphere.cpp +++ b/src/FDM/JSBSim/models/FGAtmosphere.cpp @@ -71,13 +71,13 @@ FGAtmosphere::FGAtmosphere(FGFDMExec* fdmex) : FGModel(fdmex) h = 0.0; psiw = 0.0; htab[0]=0; - htab[1]=36089.239; - htab[2]=65616.798; - htab[3]=104986.878; - htab[4]=154199.475; - htab[5]=170603.675; - htab[6]=200131.234; - htab[7]=259186.352; //ft. + htab[1]= 36089.0; + htab[2]= 65617.0; + htab[3]=104987.0; + htab[4]=154199.0; + htab[5]=167322.0; + htab[6]=232940.0; + htab[7]=278385.0; //ft. MagnitudedAccelDt = MagnitudeAccel = Magnitude = 0.0; SetTurbType( ttCulp ); @@ -86,6 +86,8 @@ FGAtmosphere::FGAtmosphere(FGFDMExec* fdmex) : FGModel(fdmex) Rhythmicity = 0.1; spike = target_time = strength = 0.0; wind_from_clockwise = 0.0; + SutherlandConstant = 198.72; // deg Rankine + Beta = 2.269690E-08; // slug/(sec ft R^0.5) T_dev_sl = T_dev = delta_T = 0.0; StandardTempOnly = false; @@ -174,52 +176,57 @@ void FGAtmosphere::Calculate(double altitude) } switch(i) { - case 1: // 36089 ft. + case 0: // Sea level + slope = -0.00356616; // R/ft. + reftemp = 518.67; // in degrees Rankine, 288.15 Kelvin + refpress = 2116.22; // psf + //refdens = 0.00237767; // slugs/cubic ft. + break; + case 1: // 36089 ft. or 11 km slope = 0; - reftemp = 389.97; - refpress = 472.452; + reftemp = 389.97; // in degrees Rankine, 216.65 Kelvin + refpress = 472.763; //refdens = 0.000706032; break; - case 2: // 65616 ft. + case 2: // 65616 ft. or 20 km slope = 0.00054864; - reftemp = 389.97; + reftemp = 389.97; // in degrees Rankine, 216.65 Kelvin refpress = 114.636; //refdens = 0.000171306; break; - case 3: // 104986 ft. - slope = 0.00153619; - reftemp = 411.57; - refpress = 8.36364; + case 3: // 104986 ft. or 32 km + slope = 0.001536192; + reftemp = 411.57; // in degrees Rankine, 228.65 Kelvin + refpress = 18.128; //refdens = 1.18422e-05; break; - case 4: // 154199 ft. + case 4: // 154199 ft. 47 km slope = 0; - reftemp = 487.17; - refpress = 0.334882; + reftemp = 487.17; // in degrees Rankine, 270.65 Kelvin + refpress = 2.316; //refdens = 4.00585e-7; break; - case 5: // 170603 ft. - slope = -0.00109728; - reftemp = 487.17; - refpress = 0.683084; + case 5: // 167322 ft. or 51 km + slope = -0.001536192; + reftemp = 487.17; // in degrees Rankine, 270.65 Kelvin + refpress = 1.398; //refdens = 8.17102e-7; break; - case 6: // 200131 ft. - slope = -0.00219456; - reftemp = 454.17; - refpress = 0.00684986; + case 6: // 232940 ft. or 71 km + slope = -0.00109728; + reftemp = 386.368; // in degrees Rankine, 214.649 Kelvin + refpress = 0.0826; //refdens = 8.77702e-9; break; - case 7: // 259186 ft. + case 7: // 278385 ft. or 84.8520 km slope = 0; - reftemp = 325.17; - refpress = 0.000122276; + reftemp = 336.5; // in degrees Rankine, 186.94 Kelvin + refpress = 0.00831; //refdens = 2.19541e-10; break; - case 0: default: // sea level slope = -0.00356616; // R/ft. - reftemp = 518.67; // R + reftemp = 518.67; // in degrees Rankine, 288.15 Kelvin refpress = 2116.22; // psf //refdens = 0.00237767; // slugs/cubic ft. break; @@ -250,7 +257,7 @@ void FGAtmosphere::Calculate(double altitude) intPressure = refpress*pow(intTemperature/reftemp,-Inertial->SLgravity()/(slope*Reng)); intDensity = intPressure/(Reng*intTemperature); } - + lastIndex=i; } @@ -263,7 +270,7 @@ void FGAtmosphere::CalculateDerived(void) T_dev = (*temperature) - GetTemperature(h); density_altitude = h + T_dev * 66.7; - if (turbType == ttStandard || ttCulp) Turbulence(); + if (turbType != ttNone) Turbulence(); vTotalWindNED = vWindNED + vGustNED + vTurbulenceNED; @@ -271,6 +278,9 @@ void FGAtmosphere::CalculateDerived(void) if (psiw < 0) psiw += 2*M_PI; soundspeed = sqrt(SHRatio*Reng*(*temperature)); + + intViscosity = Beta * pow(intTemperature, 1.5) / (SutherlandConstant + intTemperature); + intKinematicViscosity = intViscosity / intDensity; } diff --git a/src/FDM/JSBSim/models/FGAtmosphere.h b/src/FDM/JSBSim/models/FGAtmosphere.h index ed7bfcf93..b3f54f805 100644 --- a/src/FDM/JSBSim/models/FGAtmosphere.h +++ b/src/FDM/JSBSim/models/FGAtmosphere.h @@ -100,6 +100,10 @@ public: double GetDensity(double altitude); /// Returns the speed of sound in ft/sec. double GetSoundSpeed(void) const {return soundspeed;} + /// Returns the absolute viscosity. + double GetAbsoluteViscosity(void) const {return intViscosity;} + /// Returns the kinematic viscosity. + double GetKinematicViscosity(void) const {return intKinematicViscosity;} /// Returns the sea level temperature in degrees Rankine. double GetTemperatureSL(void) const { return SLtemperature; } @@ -237,6 +241,7 @@ protected: bool useExternal; double exTemperature,exDensity,exPressure; double intTemperature, intDensity, intPressure; + double SutherlandConstant, Beta, intViscosity, intKinematicViscosity; double T_dev_sl, T_dev, delta_T, density_altitude; atmType atmosphere; bool StandardTempOnly; diff --git a/src/FDM/JSBSim/models/FGAuxiliary.cpp b/src/FDM/JSBSim/models/FGAuxiliary.cpp index a04ff725d..3fa74ead7 100755 --- a/src/FDM/JSBSim/models/FGAuxiliary.cpp +++ b/src/FDM/JSBSim/models/FGAuxiliary.cpp @@ -71,6 +71,7 @@ FGAuxiliary::FGAuxiliary(FGFDMExec* fdmex) : FGModel(fdmex) qbar = 0; qbarUW = 0.0; qbarUV = 0.0; + Re = 0.0; Mach = 0.0; alpha = beta = 0.0; adot = bdot = 0.0; @@ -79,6 +80,7 @@ FGAuxiliary::FGAuxiliary(FGFDMExec* fdmex) : FGModel(fdmex) day_of_year = 1; seconds_in_day = 0.0; hoverbmac = hoverbcg = 0.0; + tatc = RankineToCelsius(tat); vPilotAccel.InitMatrix(); vPilotAccelN.InitMatrix(); @@ -200,6 +202,8 @@ bool FGAuxiliary::Run() alpha = beta = adot = bdot = 0; } + Re = Vt * Aircraft->Getcbar() / Atmosphere->GetKinematicViscosity(); + qbar = 0.5*Atmosphere->GetDensity()*Vt*Vt; qbarUW = 0.5*Atmosphere->GetDensity()*(vAeroUVW(eU)*vAeroUVW(eU) + vAeroUVW(eW)*vAeroUVW(eW)); qbarUV = 0.5*Atmosphere->GetDensity()*(vAeroUVW(eU)*vAeroUVW(eU) + vAeroUVW(eV)*vAeroUVW(eV)); @@ -349,6 +353,7 @@ void FGAuxiliary::bind(void) PropertyManager->Tie("aero/alpha-deg", this, inDegrees, (PMF)&FGAuxiliary::Getalpha); PropertyManager->Tie("aero/beta-deg", this, inDegrees, (PMF)&FGAuxiliary::Getbeta); PropertyManager->Tie("aero/mag-beta-deg", this, inDegrees, (PMF)&FGAuxiliary::GetMagBeta); + PropertyManager->Tie("aero/Re", this, &FGAuxiliary::GetReynoldsNumber); PropertyManager->Tie("aero/qbar-psf", this, &FGAuxiliary::Getqbar, &FGAuxiliary::Setqbar, true); PropertyManager->Tie("aero/qbarUW-psf", this, &FGAuxiliary::GetqbarUW, &FGAuxiliary::SetqbarUW, true); PropertyManager->Tie("aero/qbarUV-psf", this, &FGAuxiliary::GetqbarUV, &FGAuxiliary::SetqbarUV, true); diff --git a/src/FDM/JSBSim/models/FGAuxiliary.h b/src/FDM/JSBSim/models/FGAuxiliary.h index e551f5e3d..8d47ce8cc 100644 --- a/src/FDM/JSBSim/models/FGAuxiliary.h +++ b/src/FDM/JSBSim/models/FGAuxiliary.h @@ -181,14 +181,15 @@ public: double GetMagBeta (int unit) const { if (unit == inDegrees) return fabs(beta)*radtodeg; else cerr << "Bad units" << endl; return 0.0;} - double Getqbar (void) const { return qbar; } - double GetqbarUW (void) const { return qbarUW; } - double GetqbarUV (void) const { return qbarUV; } - double GetVt (void) const { return Vt; } - double GetVground (void) const { return Vground; } - double GetMach (void) const { return Mach; } - double GetMachU (void) const { return MachU; } - double GetNz (void) const { return Nz; } + double Getqbar (void) const { return qbar; } + double GetqbarUW (void) const { return qbarUW; } + double GetqbarUV (void) const { return qbarUV; } + double GetReynoldsNumber(void) const { return Re; } + double GetVt (void) const { return Vt; } + double GetVground (void) const { return Vground; } + double GetMach (void) const { return Mach; } + double GetMachU (void) const { return MachU; } + double GetNz (void) const { return Nz; } double GetHOverBCG(void) const { return hoverbcg; } double GetHOverBMAC(void) const { return hoverbmac; } @@ -247,6 +248,7 @@ private: double Vt, Vground, Mach, MachU; double qbar, qbarUW, qbarUV; + double Re; // Reynolds Number = V*c/mu double alpha, beta; double adot,bdot; double psigt, gamma; diff --git a/src/FDM/JSBSim/models/FGFCS.cpp b/src/FDM/JSBSim/models/FGFCS.cpp index 95abe9ebb..3d974bf6e 100644 --- a/src/FDM/JSBSim/models/FGFCS.cpp +++ b/src/FDM/JSBSim/models/FGFCS.cpp @@ -46,7 +46,6 @@ INCLUDES #include #include #include -#include #include #include #include @@ -514,6 +513,10 @@ bool FGFCS::Load(Element* el, SystemType systype) return false; } else { document = LoadXMLDocument(file); + if (!document) { + cerr << "Error loading file " << file << endl; + return false; + } name = document->GetAttributeValue("name"); } } else { @@ -551,7 +554,7 @@ bool FGFCS::Load(Element* el, SystemType systype) } // After reading interface properties in a file, read properties in the local - // flight_control, autopiot, or system element. This allows general-purpose + // flight_control, autopilot, or system element. This allows general-purpose // systems to be defined in a file, with overrides or initial loaded constants // supplied in the relevant element of the aircraft configuration file. @@ -967,7 +970,7 @@ void FGFCS::Debug(int from) if (debug_lvl & 1) { // Standard console startup message output if (from == 2) { // Loader - cout << endl << " Flight Control (" << Name << ")" << endl; + cout << endl << " " << Name << endl; } } if (debug_lvl & 2 ) { // Instantiation/Destruction notification diff --git a/src/FDM/JSBSim/models/FGGroundReactions.cpp b/src/FDM/JSBSim/models/FGGroundReactions.cpp index f8067cfba..4638c1f39 100644 --- a/src/FDM/JSBSim/models/FGGroundReactions.cpp +++ b/src/FDM/JSBSim/models/FGGroundReactions.cpp @@ -89,17 +89,14 @@ bool FGGroundReactions::Run(void) vForces.InitMatrix(); vMoments.InitMatrix(); - // Sum forces and moments for all gear, here. - // Some optimizations may be made here - or rather in the gear code itself. - // The gear ::Run() method is called several times - once for each gear. - // Perhaps there is some commonality for things which only need to be - // calculated once. - if ( Propagate->GetDistanceAGL() < 300.0 ) { // Only execute gear code below 300 feet - for (unsigned int i=0; iForce(); - vMoments += lGear[i]->Moment(); - } - + // Sum forces and moments for all gear, here. + // Some optimizations may be made here - or rather in the gear code itself. + // The gear ::Run() method is called several times - once for each gear. + // Perhaps there is some commonality for things which only need to be + // calculated once. + for (unsigned int i=0; iForce(); + vMoments += lGear[i]->Moment(); } return false; diff --git a/src/FDM/JSBSim/models/FGInertial.cpp b/src/FDM/JSBSim/models/FGInertial.cpp index 75a699f99..39aa91ddc 100644 --- a/src/FDM/JSBSim/models/FGInertial.cpp +++ b/src/FDM/JSBSim/models/FGInertial.cpp @@ -54,7 +54,7 @@ FGInertial::FGInertial(FGFDMExec* fgex) : FGModel(fgex) { Name = "FGInertial"; - // Defaults + // Earth defaults RotationRate = 0.00007292115; GM = 14.07644180E15; // WGS84 value RadiusReference = 20925650.00; // Equatorial radius (WGS84) @@ -64,6 +64,18 @@ FGInertial::FGInertial(FGFDMExec* fgex) : FGModel(fgex) b = 20855486.5951; // WGS84 semiminor axis length in feet earthPosAngle = 0.0; + // Lunar defaults + /* + RotationRate = 0.0000026617; + GM = 1.7314079E14; // Lunar GM + 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 + b = 5695439.63; // semiminor axis length in feet + earthPosAngle = 0.0; + */ + gAccelReference = GM/(RadiusReference*RadiusReference); gAccel = GM/(RadiusReference*RadiusReference); diff --git a/src/FDM/JSBSim/models/FGLGear.cpp b/src/FDM/JSBSim/models/FGLGear.cpp index c1110dca6..080111976 100644 --- a/src/FDM/JSBSim/models/FGLGear.cpp +++ b/src/FDM/JSBSim/models/FGLGear.cpp @@ -57,8 +57,9 @@ static const char *IdHdr = ID_LGEAR; CLASS IMPLEMENTATION %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/ -FGLGear::FGLGear(Element* el, FGFDMExec* fdmex, int number) : Exec(fdmex), - GearNumber(number) +FGLGear::FGLGear(Element* el, FGFDMExec* fdmex, int number) : + GearNumber(number), + Exec(fdmex) { Element *force_table=0; Element *dampCoeff=0; @@ -341,97 +342,96 @@ FGColumnVector3& FGLGear::Force(void) if (isRetractable) ComputeRetractionState(); - if (!GearDown) return vForce; // return the null vForce column vector + if (GearDown) { - vWhlBodyVec = MassBalance->StructuralToBody(vXYZ); // Get wheel in body frame - vLocalGear = Propagate->GetTb2l() * vWhlBodyVec; // Get local frame wheel location + vWhlBodyVec = MassBalance->StructuralToBody(vXYZ); // Get wheel in body frame + vLocalGear = Propagate->GetTb2l() * vWhlBodyVec; // Get local frame wheel location - gearLoc = Propagate->GetLocation().LocalToLocation(vLocalGear); - compressLength = -Exec->GetGroundCallback()->GetAGLevel(t, gearLoc, contact, normal, cvel); + gearLoc = Propagate->GetLocation().LocalToLocation(vLocalGear); + compressLength = -Exec->GetGroundCallback()->GetAGLevel(t, gearLoc, contact, normal, cvel); - // The compression length is measured in the Z-axis, only, at this time. + // The compression length is measured in the Z-axis, only, at this time. - if (compressLength > 0.00) { + if (compressLength > 0.00) { - WOW = true; + WOW = true; - // [The next equation should really use the vector to the contact patch of - // the tire including the strut compression and not the original vWhlBodyVec.] + // [The next equation should really use the vector to the contact patch of + // the tire including the strut compression and not the original vWhlBodyVec.] - vWhlVelVec = Propagate->GetTb2l() * (Propagate->GetPQR() * vWhlBodyVec); - vWhlVelVec += Propagate->GetVel() - cvel; - compressSpeed = vWhlVelVec(eZ); + vWhlVelVec = Propagate->GetTb2l() * (Propagate->GetPQR() * vWhlBodyVec); + vWhlVelVec += Propagate->GetVel() - cvel; + compressSpeed = vWhlVelVec(eZ); - InitializeReporting(); - ComputeBrakeForceCoefficient(); - ComputeSteeringAngle(); - ComputeSlipAngle(); - ComputeSideForceCoefficient(); - ComputeVerticalStrutForce(); + InitializeReporting(); + ComputeBrakeForceCoefficient(); + ComputeSteeringAngle(); + ComputeSlipAngle(); + ComputeSideForceCoefficient(); + ComputeVerticalStrutForce(); - // Compute the forces in the wheel ground plane. + // Compute the forces in the wheel ground plane. - double sign = RollingWhlVel>0?1.0:(RollingWhlVel<0?-1.0:0.0); - RollingForce = ((1.0 - TirePressureNorm) * 30 + vLocalForce(eZ) * BrakeFCoeff) * sign; - SideForce = vLocalForce(eZ) * FCoeff; + double sign = RollingWhlVel>0?1.0:(RollingWhlVel<0?-1.0:0.0); + RollingForce = ((1.0 - TirePressureNorm) * 30 + vLocalForce(eZ) * BrakeFCoeff) * sign; + SideForce = vLocalForce(eZ) * FCoeff; - // Transform these forces back to the local reference frame. + // Transform these forces back to the local reference frame. - vLocalForce(eX) = RollingForce*CosWheel - SideForce*SinWheel; - vLocalForce(eY) = SideForce*CosWheel + RollingForce*SinWheel; + vLocalForce(eX) = RollingForce*CosWheel - SideForce*SinWheel; + vLocalForce(eY) = SideForce*CosWheel + RollingForce*SinWheel; - // Transform the forces back to the body frame and compute the moment. + // Transform the forces back to the body frame and compute the moment. - vForce = Propagate->GetTl2b() * vLocalForce; + vForce = Propagate->GetTl2b() * vLocalForce; -// Start experimental section for gear jitter reduction -// -// Lag and attenuate the XY-plane forces dependent on velocity + // Lag and attenuate the XY-plane forces dependent on velocity - double ca, cb, denom; - FGColumnVector3 Output; + double ca, cb, denom; + FGColumnVector3 Output; -// This code implements a lag filter, C/(s + C) where -// "C" is the filter coefficient. When "C" is chosen at the -// frame rate (in Hz), the jittering is significantly reduced. This is because -// the jitter is present *at* the execution rate. -// If a coefficient is set to something equal to or less than zero, the filter -// is bypassed. + // This code implements a lag filter, C/(s + C) where + // "C" is the filter coefficient. When "C" is chosen at the + // frame rate (in Hz), the jittering is significantly reduced. This is because + // the jitter is present *at* the execution rate. + // If a coefficient is set to something equal to or less than zero, the filter + // is bypassed. - if (LongForceLagFilterCoeff > 0) { - denom = 2.00 + dT*LongForceLagFilterCoeff; - ca = dT*LongForceLagFilterCoeff / denom; - cb = (2.00 - dT*LongForceLagFilterCoeff) / denom; - Output(eX) = vForce(eX) * ca + prevIn(eX) * ca + prevOut(eX) * cb; - vForce(eX) = Output(eX); + if (LongForceLagFilterCoeff > 0) { + denom = 2.00 + dT*LongForceLagFilterCoeff; + ca = dT*LongForceLagFilterCoeff / denom; + cb = (2.00 - dT*LongForceLagFilterCoeff) / denom; + Output(eX) = vForce(eX) * ca + prevIn(eX) * ca + prevOut(eX) * cb; + vForce(eX) = Output(eX); + } + if (LatForceLagFilterCoeff > 0) { + denom = 2.00 + dT*LatForceLagFilterCoeff; + ca = dT*LatForceLagFilterCoeff / denom; + cb = (2.00 - dT*LatForceLagFilterCoeff) / denom; + Output(eY) = vForce(eY) * ca + prevIn(eY) * ca + prevOut(eY) * cb; + vForce(eY) = Output(eY); + } + + prevIn = vForce; + prevOut = Output; + + if ((fabs(RollingWhlVel) <= RFRV) && RFRV > 0) vForce(eX) *= fabs(RollingWhlVel)/RFRV; + if ((fabs(SideWhlVel) <= SFRV) && SFRV > 0) vForce(eY) *= fabs(SideWhlVel)/SFRV; + + // End section for attentuating gear jitter + + vMoment = vWhlBodyVec * vForce; + + } else { // Gear is NOT compressed + + WOW = false; + compressLength = 0.0; + + // Return to neutral position between 1.0 and 0.8 gear pos. + SteerAngle *= max(GetGearUnitPos()-0.8, 0.0)/0.2; + + ResetReporting(); } - if (LatForceLagFilterCoeff > 0) { - denom = 2.00 + dT*LatForceLagFilterCoeff; - ca = dT*LatForceLagFilterCoeff / denom; - cb = (2.00 - dT*LatForceLagFilterCoeff) / denom; - Output(eY) = vForce(eY) * ca + prevIn(eY) * ca + prevOut(eY) * cb; - vForce(eY) = Output(eY); - } - - prevIn = vForce; - prevOut = Output; - - if ((fabs(RollingWhlVel) <= RFRV) && RFRV > 0) vForce(eX) *= fabs(RollingWhlVel)/RFRV; - if ((fabs(SideWhlVel) <= SFRV) && SFRV > 0) vForce(eY) *= fabs(SideWhlVel)/SFRV; - -// End section for attentuating gear jitter - - vMoment = vWhlBodyVec * vForce; - - } else { // Gear is NOT compressed - - WOW = false; - compressLength = 0.0; - - // Return to neutral position between 1.0 and 0.8 gear pos. - SteerAngle *= max(GetGearUnitPos()-0.8, 0.0)/0.2; - - ResetReporting(); } ReportTakeoffOrLanding(); @@ -452,6 +452,7 @@ void FGLGear::ComputeRetractionState(void) double gearPos = GetGearUnitPos(); if (gearPos < 0.01) { GearUp = true; + WOW = false; GearDown = false; } else if (gearPos > 0.99) { GearDown = true; @@ -575,7 +576,7 @@ void FGLGear::ReportTakeoffOrLanding(void) LandingDistanceTraveled += Auxiliary->GetVground()*deltaT; if (StartedGroundRun) { - TakeoffDistanceTraveled50ft += Auxiliary->GetVground()*deltaT; + TakeoffDistanceTraveled50ft += Auxiliary->GetVground()*deltaT; if (WOW) TakeoffDistanceTraveled += Auxiliary->GetVground()*deltaT; } @@ -739,6 +740,8 @@ void FGLGear::bind(void) snprintf(property_name, 80, "gear/unit[%d]/z-position", GearNumber); Exec->GetPropertyManager()->Tie( property_name, (FGLGear*)this, &FGLGear::GetZPosition, &FGLGear::SetZPosition); + snprintf(property_name, 80, "gear/unit[%d]/compression-ft", GearNumber); + Exec->GetPropertyManager()->Tie( property_name, &compressLength ); } if( isRetractable ) { @@ -752,6 +755,8 @@ void FGLGear::bind(void) void FGLGear::Report(ReportType repType) { + if (fabs(TakeoffDistanceTraveled) < 0.001) return; // Don't print superfluous reports + switch(repType) { case erLand: cout << endl << "Touchdown report for " << name << endl; diff --git a/src/FDM/JSBSim/models/FGLGear.h b/src/FDM/JSBSim/models/FGLGear.h index 3dc81b7c1..4de980edf 100644 --- a/src/FDM/JSBSim/models/FGLGear.h +++ b/src/FDM/JSBSim/models/FGLGear.h @@ -257,6 +257,9 @@ public: /// Sets the brake value in percent (0 - 100) inline void SetBrake(double bp) {brakePct = bp;} + /// Sets the weight-on-wheels flag. + void SetWOW(bool wow) {WOW = wow;} + /** Set the console touchdown reporting feature @param flag true turns on touchdown reporting, false turns it off */ inline void SetReport(bool flag) { ReportEnable = flag; } diff --git a/src/FDM/JSBSim/models/FGOutput.cpp b/src/FDM/JSBSim/models/FGOutput.cpp index 032245e85..4b3a95f88 100644 --- a/src/FDM/JSBSim/models/FGOutput.cpp +++ b/src/FDM/JSBSim/models/FGOutput.cpp @@ -247,6 +247,7 @@ void FGOutput::DelimitedOutput(string fname) if (SubSystems & ssVelocities) { outstream << delimeter; outstream << "q bar (psf)" + delimeter; + outstream << "Reynolds Number" + delimeter; outstream << "V_{Total} (ft/s)" + delimeter; outstream << "V_{Inertial} (ft/s)" + delimeter; outstream << "UBody" + delimeter + "VBody" + delimeter + "WBody" + delimeter; @@ -266,6 +267,9 @@ void FGOutput::DelimitedOutput(string fname) if (SubSystems & ssAtmosphere) { outstream << delimeter; outstream << "Rho (slugs/ft^3)" + delimeter; + outstream << "Absolute Viscosity" + delimeter; + outstream << "Kinematic Viscosity" + delimeter; + outstream << "Temperature (R)" + delimeter; outstream << "P_{SL} (psf)" + delimeter; outstream << "P_{Ambient} (psf)" + delimeter; outstream << "Turbulence Magnitude (ft/sec)" + delimeter; @@ -274,17 +278,17 @@ void FGOutput::DelimitedOutput(string fname) } if (SubSystems & ssMassProps) { outstream << delimeter; - outstream << "I_xx" + delimeter; - outstream << "I_xy" + delimeter; - outstream << "I_xz" + delimeter; - outstream << "I_yx" + delimeter; - outstream << "I_yy" + delimeter; - outstream << "I_yz" + delimeter; - outstream << "I_zx" + delimeter; - outstream << "I_zy" + delimeter; - outstream << "I_zz" + delimeter; + outstream << "I_{xx}" + delimeter; + outstream << "I_{xy}" + delimeter; + outstream << "I_{xz}" + delimeter; + outstream << "I_{yx}" + delimeter; + outstream << "I_{yy}" + delimeter; + outstream << "I_{yz}" + delimeter; + outstream << "I_{zx}" + delimeter; + outstream << "I_{zy}" + delimeter; + outstream << "I_{zz}" + delimeter; outstream << "Mass" + delimeter; - outstream << "X_cg" + delimeter + "Y_cg" + delimeter + "Z_cg"; + outstream << "X_{cg}" + delimeter + "Y_{cg}" + delimeter + "Z_{cg}"; } if (SubSystems & ssPropagate) { outstream << delimeter; @@ -348,6 +352,7 @@ void FGOutput::DelimitedOutput(string fname) if (SubSystems & ssVelocities) { outstream << delimeter; outstream << Auxiliary->Getqbar() << delimeter; + outstream << Auxiliary->GetReynoldsNumber() << delimeter; outstream << setprecision(12) << Auxiliary->GetVt() << delimeter; outstream << Propagate->GetInertialVelocityMagnitude() << delimeter; outstream << setprecision(12) << Propagate->GetUVW().Dump(delimeter) << delimeter; @@ -367,6 +372,9 @@ void FGOutput::DelimitedOutput(string fname) if (SubSystems & ssAtmosphere) { outstream << delimeter; outstream << Atmosphere->GetDensity() << delimeter; + outstream << Atmosphere->GetAbsoluteViscosity() << delimeter; + outstream << Atmosphere->GetKinematicViscosity() << delimeter; + outstream << Atmosphere->GetTemperature() << delimeter; outstream << Atmosphere->GetPressureSL() << delimeter; outstream << Atmosphere->GetPressure() << delimeter; outstream << Atmosphere->GetTurbMagnitude() << delimeter; diff --git a/src/FDM/JSBSim/models/FGPropulsion.cpp b/src/FDM/JSBSim/models/FGPropulsion.cpp index 2ec2e2bb3..82092f8ad 100644 --- a/src/FDM/JSBSim/models/FGPropulsion.cpp +++ b/src/FDM/JSBSim/models/FGPropulsion.cpp @@ -184,7 +184,7 @@ bool FGPropulsion::GetSteadyState(void) while (!steady && j < 6000) { Engines[i]->Calculate(); lastThrust = currentThrust; - currentThrust = Engines[i]->GetThrust(); + currentThrust = Engines[i]->GetThruster()->GetThrust(); if (fabs(lastThrust-currentThrust) < 0.0001) { steady_count++; if (steady_count > 120) { diff --git a/src/FDM/JSBSim/models/atmosphere/FGMSIS.cpp b/src/FDM/JSBSim/models/atmosphere/FGMSIS.cpp index 6eeff5af8..e3197404d 100755 --- a/src/FDM/JSBSim/models/atmosphere/FGMSIS.cpp +++ b/src/FDM/JSBSim/models/atmosphere/FGMSIS.cpp @@ -136,7 +136,7 @@ bool MSIS::InitModel(void) pressure = &intPressure; density = &intDensity; - useExternal=false; + UseInternal(); return true; } @@ -174,19 +174,11 @@ bool MSIS::Run(void) intTemperature = output.t[1] * 1.8; intDensity = output.d[5] * 1.940321; intPressure = 1716.488 * intDensity * intTemperature; - soundspeed = sqrt(2403.0832 * intTemperature); //cout << "T=" << intTemperature << " D=" << intDensity << " P="; //cout << intPressure << " a=" << soundspeed << endl; } - if (turbType != ttNone) { - Turbulence(); - vWindNED += vTurbulenceNED; - } - - if (vWindNED(1) != 0.0) psiw = atan2( vWindNED(2), vWindNED(1) ); - - if (psiw < 0) psiw += 2*M_PI; + CalculateDerived(); Debug(2); diff --git a/src/FDM/JSBSim/models/propulsion/FGEngine.cpp b/src/FDM/JSBSim/models/propulsion/FGEngine.cpp index 49f67e2f4..49937ab8f 100644 --- a/src/FDM/JSBSim/models/propulsion/FGEngine.cpp +++ b/src/FDM/JSBSim/models/propulsion/FGEngine.cpp @@ -119,7 +119,7 @@ FGEngine::FGEngine(FGFDMExec* exec, Element* engine_element, int engine_number) snprintf(property_name, 80, "propulsion/engine[%d]/set-running", EngineNumber); PropertyManager->Tie( property_name, this, &FGEngine::GetRunning, &FGEngine::SetRunning ); snprintf(property_name, 80, "propulsion/engine[%u]/thrust-lbs", EngineNumber); - PropertyManager->Tie( property_name, this, &FGEngine::GetThrust); + PropertyManager->Tie( property_name, Thruster, &FGThruster::GetThrust); snprintf(property_name, 80, "propulsion/engine[%u]/fuel-flow-rate-pps", EngineNumber); PropertyManager->Tie( property_name, this, &FGEngine::GetFuelFlowRate); @@ -138,7 +138,6 @@ FGEngine::~FGEngine() void FGEngine::ResetToIC(void) { - Thrust = 0.0; Throttle = 0.0; Mixture = 1.0; Starter = false; diff --git a/src/FDM/JSBSim/models/propulsion/FGEngine.h b/src/FDM/JSBSim/models/propulsion/FGEngine.h index dd8b1302b..783fbc0fe 100644 --- a/src/FDM/JSBSim/models/propulsion/FGEngine.h +++ b/src/FDM/JSBSim/models/propulsion/FGEngine.h @@ -150,7 +150,6 @@ public: virtual double getFuelFlow_gph () const {return FuelFlow_gph;} virtual double getFuelFlow_pph () const {return FuelFlow_pph;} virtual double GetFuelFlowRate(void) const {return FuelFlowRate;} - virtual double GetThrust(void) const { return Thrust; } virtual bool GetStarved(void) { return Starved; } virtual bool GetRunning(void) const { return Running; } virtual bool GetCranking(void) { return Cranking; } @@ -216,7 +215,6 @@ protected: double MaxThrottle; double MinThrottle; - double Thrust; double Throttle; double Mixture; double FuelExpended; diff --git a/src/FDM/JSBSim/models/propulsion/FGPiston.cpp b/src/FDM/JSBSim/models/propulsion/FGPiston.cpp index e23d981a8..bac4eaa91 100644 --- a/src/FDM/JSBSim/models/propulsion/FGPiston.cpp +++ b/src/FDM/JSBSim/models/propulsion/FGPiston.cpp @@ -132,21 +132,6 @@ FGPiston::FGPiston(FGFDMExec* exec, Element* el, int engine_number) *Lookup_Combustion_Efficiency << 1.60 << 0.525; *Lookup_Combustion_Efficiency << 2.00 << 0.345; - Power_Mixture_Correlation = new FGTable(13); - *Power_Mixture_Correlation << (14.7/1.6) << 0.780; - *Power_Mixture_Correlation << 10 << 0.860; - *Power_Mixture_Correlation << 11 << 0.935; - *Power_Mixture_Correlation << 12 << 0.980; - *Power_Mixture_Correlation << 13 << 1.000; - *Power_Mixture_Correlation << 14 << 0.990; - *Power_Mixture_Correlation << 15 << 0.964; - *Power_Mixture_Correlation << 16 << 0.925; - *Power_Mixture_Correlation << 17 << 0.880; - *Power_Mixture_Correlation << 18 << 0.830; - *Power_Mixture_Correlation << 19 << 0.785; - *Power_Mixture_Correlation << 20 << 0.740; - *Power_Mixture_Correlation << (14.7/0.6) << 0.58; - Mixture_Efficiency_Correlation = new FGTable(15); *Mixture_Efficiency_Correlation << 0.05000 << 0.00000; *Mixture_Efficiency_Correlation << 0.05137 << 0.00862; @@ -165,21 +150,6 @@ FGPiston::FGPiston(FGFDMExec* exec, Element* el, int engine_number) *Mixture_Efficiency_Correlation << 0.12500 << 0.00000; -/* -Manifold_Pressure_Lookup = new - - 0 0.2 0.4 0.6 0.8 1 -0 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 -1000 0.7778 0.8212 0.8647 0.9081 0.9516 0.9950 -2000 0.5556 0.6424 0.7293 0.8162 0.9031 0.9900 -3000 0.3333 0.4637 0.5940 0.7243 0.8547 0.9850 -4000 0.2000 0.2849 0.4587 0.6324 0.8062 0.9800 -5000 0.2000 0.2000 0.3233 0.5406 0.7578 0.9750 -6000 0.2000 0.2000 0.2000 0.4487 0.7093 0.9700 -7000 0.2000 0.2000 0.2000 0.2000 0.4570 0.7611 -8000 0.2000 0.2000 0.2000 0.2000 0.2047 0.5522 -*/ - // Read inputs from engine data file where present. if (el->FindElement("minmp")) // Should have ELSE statement telling default value used? @@ -238,10 +208,17 @@ Manifold_Pressure_Lookup = new RatedAltitude[2] = el->FindElementValueAsNumberConvertTo("ratedaltitude3", "FT"); } + MaxManifoldPressure_Percent = MaxManifoldPressure_inHg / 29.92; // Create a BSFC to match the engine if not provided if (BSFC < 0) { BSFC = ( Displacement * MaxRPM * volumetric_efficiency ) / (9411 * MaxHP); + BSFC *= (MaxManifoldPressure_Percent * MaxManifoldPressure_Percent * MaxManifoldPressure_Percent); } + if ( MaxManifoldPressure_inHg > 29.9 ) { // Don't allow boosting with a bogus number + MaxManifoldPressure_inHg = 29.9; + MaxManifoldPressure_Percent = MaxManifoldPressure_inHg / 29.92; + } + char property_name[80]; snprintf(property_name, 80, "propulsion/engine[%d]/power-hp", EngineNumber); PropertyManager->Tie(property_name, &HP); @@ -299,7 +276,9 @@ Manifold_Pressure_Lookup = new BoostSpeed = 0; } bBoostOverride = (BoostOverride == 1 ? true : false); - if (MinThrottle < 0.001) MinThrottle = 0.001; //MinThrottle is a denominator in a power equation so it can't be zero + if (MinThrottle < 0.12) MinThrottle = 0.12; //MinThrottle is limited to 0.12 to prevent the + // throttle area equation from going negative + // 0.12 is 1% of maximum area Debug(0); // Call Debug() routine from constructor if needed } @@ -308,7 +287,6 @@ Manifold_Pressure_Lookup = new FGPiston::~FGPiston() { delete Lookup_Combustion_Efficiency; - delete Power_Mixture_Correlation; delete Mixture_Efficiency_Correlation; Debug(1); // Call Debug() routine from constructor if needed } @@ -338,7 +316,8 @@ double FGPiston::Calculate(void) if (FuelFlow_gph > 0.0) ConsumeFuel(); Throttle = FCS->GetThrottlePos(EngineNumber); - ThrottlePos = MinThrottle+((MaxThrottle-MinThrottle)*Throttle ); + // calculate the throttle plate angle. 1 unit is pi/2 radians. + ThrottleAngle = MinThrottle+((MaxThrottle-MinThrottle)*Throttle ); Mixture = FCS->GetMixturePos(EngineNumber); // @@ -367,8 +346,7 @@ double FGPiston::Calculate(void) // Running = false; doEnginePower(); -if(HP<0.1250) - Running = false; + if (HP < 0.1250) Running = false; doEGT(); doCHT(); @@ -390,8 +368,6 @@ if(HP<0.1250) double FGPiston::CalcFuelNeed(void) { double dT = State->Getdt() * Propulsion->GetRate(); - FuelFlow_pph = FuelFlow_gph * 6.0; // Assumes 6 lbs / gallon - FuelFlowRate = FuelFlow_pph / 3600.0; FuelExpended = FuelFlowRate * dT; return FuelExpended; } @@ -510,20 +486,18 @@ void FGPiston::doBoostControl(void) * from the throttle position, turbo/supercharger boost control * system, engine speed and local ambient air density. * - * TODO: changes in MP should not be instantaneous -- introduce - * a lag between throttle changes and MP changes, to allow pressure - * to build up or disperse. - * - * Inputs: minMAP, maxMAP, p_amb, Throttle + * Inputs: p_amb, Throttle, MaxManifoldPressure_Percent, ThrottleAngle + * RPM, MaxRPM * * Outputs: MAP, ManifoldPressure_inHg */ void FGPiston::doMAP(void) { - suction_loss = RPM > 0.0 ? ThrottlePos * MaxRPM / RPM : 1.0; - if (suction_loss > 1.0) suction_loss = 1.0; - MAP = p_amb * suction_loss; + // estimate throttle plate area. This maps 0.2 -> 0.1 for historical performance reasons + double throttle_area = ThrottleAngle * 1.125 - 0.125; + map_coefficient = pow ((throttle_area * MaxManifoldPressure_Percent),RPM/MaxRPM); + MAP = p_amb * map_coefficient; if(Boosted) { // If takeoff boost is fitted, we currently assume the following throttle map: @@ -549,7 +523,7 @@ void FGPiston::doMAP(void) } } // Boost the manifold pressure. - double boost_factor = BoostMul[BoostSpeed] * suction_loss * RPM/RatedRPM[BoostSpeed]; + double boost_factor = BoostMul[BoostSpeed] * map_coefficient * RPM/RatedRPM[BoostSpeed]; if (boost_factor < 1.0) boost_factor = 1.0; // boost will never reduce the MAP MAP *= boost_factor; // Now clip the manifold pressure to BCV or Wastegate setting. @@ -575,7 +549,7 @@ void FGPiston::doMAP(void) * (used in CHT calculation for air-cooled engines). * * Inputs: p_amb, R_air, T_amb, MAP, Displacement, - * RPM, volumetric_efficiency, ThrottlePos + * RPM, volumetric_efficiency, ThrottleAngle * * TODO: Model inlet manifold air temperature. * @@ -587,7 +561,7 @@ void FGPiston::doAirFlow(void) rho_air = p_amb / (R_air * T_amb); double displacement_SI = Displacement * in3tom3; double swept_volume = (displacement_SI * (RPM/60)) / 2; - double v_dot_air = swept_volume * volumetric_efficiency * suction_loss; + double v_dot_air = swept_volume * volumetric_efficiency * map_coefficient; double rho_air_manifold = MAP / (R_air * T_amb); m_dot_air = v_dot_air * rho_air_manifold; @@ -609,11 +583,9 @@ void FGPiston::doFuelFlow(void) // double AFR = 10+(12*(1-Mixture));// mixture 10:1 to 22:1 // m_dot_fuel = m_dot_air / AFR; m_dot_fuel = (m_dot_air * equivalence_ratio) / 14.7; - FuelFlow_gph = m_dot_fuel - * 3600 // seconds to hours - * 2.2046 // kg to lb - / 6.0; // lb to gal_us of gasoline -// / 6.6; // lb to gal_us of kerosene + FuelFlowRate = m_dot_fuel * 2.2046; // kg to lb + FuelFlow_pph = FuelFlowRate * 3600; // seconds to hours + FuelFlow_gph = FuelFlow_pph / 6.0; // Assumes 6 lbs / gallon } //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% @@ -649,7 +621,7 @@ void FGPiston::doEnginePower(void) if ( Magnetos != 3 ) power *= SparkFailDrop; - HP = (FuelFlow_gph * 6.0 / BSFC )* ME * suction_loss * power; + HP = (FuelFlow_gph * 6.0 / BSFC )* ME * map_coefficient * power; } else { @@ -871,17 +843,13 @@ void FGPiston::Debug(int from) cout << " IdleRPM: " << IdleRPM << endl; cout << " MaxThrottle: " << MaxThrottle << endl; cout << " MinThrottle: " << MinThrottle << endl; + cout << " BSFC: " << BSFC << endl; cout << endl; cout << " Combustion Efficiency table:" << endl; Lookup_Combustion_Efficiency->Print(); cout << endl; - cout << endl; - cout << " Power Mixture Correlation table:" << endl; - Power_Mixture_Correlation->Print(); - cout << endl; - cout << endl; cout << " Mixture Efficiency Correlation table:" << endl; Mixture_Efficiency_Correlation->Print(); diff --git a/src/FDM/JSBSim/models/propulsion/FGPiston.h b/src/FDM/JSBSim/models/propulsion/FGPiston.h index d4fa080f3..8878c5040 100644 --- a/src/FDM/JSBSim/models/propulsion/FGPiston.h +++ b/src/FDM/JSBSim/models/propulsion/FGPiston.h @@ -205,8 +205,7 @@ public: double getRPM(void) {return RPM;} protected: - double ThrottlePos; - + double ThrottleAngle; private: int crank_counter; @@ -252,6 +251,7 @@ private: // double MinManifoldPressure_inHg; // Inches Hg double MaxManifoldPressure_inHg; // Inches Hg + double MaxManifoldPressure_Percent; // MaxManifoldPressure / 29.92 double Displacement; // cubic inches double MaxHP; // horsepower double SparkFailDrop; // drop of power due to spark failure @@ -304,7 +304,7 @@ private: // double rho_air; double volumetric_efficiency; - double suction_loss; + double map_coefficient; double m_dot_air; double equivalence_ratio; double m_dot_fuel; diff --git a/src/FDM/JSBSim/models/propulsion/FGRocket.cpp b/src/FDM/JSBSim/models/propulsion/FGRocket.cpp index 66390010c..b4edabfdd 100644 --- a/src/FDM/JSBSim/models/propulsion/FGRocket.cpp +++ b/src/FDM/JSBSim/models/propulsion/FGRocket.cpp @@ -106,6 +106,7 @@ FGRocket::~FGRocket(void) double FGRocket::Calculate(void) { double dT = State->Getdt()*Propulsion->GetRate(); + double thrust; if (!Flameout && !Starved) ConsumeFuel(); @@ -135,7 +136,7 @@ double FGRocket::Calculate(void) if (Throttle < MinThrottle || Starved) { // Combustion not supported - PctPower = Thrust = 0.0; // desired thrust + PctPower = 0.0; // desired thrust Flameout = true; VacThrust = 0.0; @@ -149,10 +150,10 @@ double FGRocket::Calculate(void) } // End thrust calculations - Thrust = Thruster->Calculate(VacThrust); - It += Thrust * dT; + thrust = Thruster->Calculate(VacThrust); + It += thrust * dT; - return Thrust; + return thrust; } //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% diff --git a/src/FDM/JSBSim/models/propulsion/FGTank.cpp b/src/FDM/JSBSim/models/propulsion/FGTank.cpp index caa8ef2f6..af174d46d 100644 --- a/src/FDM/JSBSim/models/propulsion/FGTank.cpp +++ b/src/FDM/JSBSim/models/propulsion/FGTank.cpp @@ -148,6 +148,8 @@ FGTank::FGTank(FGFDMExec* exec, Element* el, int tank_number) if (Temperature != -9999.0) InitialTemperature = Temperature = FahrenheitToCelsius(Temperature); Area = 40.0 * pow(Capacity/1975, 0.666666667); + CalculateInertias(); + Debug(0); } diff --git a/src/FDM/JSBSim/models/propulsion/FGThruster.cpp b/src/FDM/JSBSim/models/propulsion/FGThruster.cpp index b802d6e8c..9bbe8ed1c 100644 --- a/src/FDM/JSBSim/models/propulsion/FGThruster.cpp +++ b/src/FDM/JSBSim/models/propulsion/FGThruster.cpp @@ -86,7 +86,7 @@ FGThruster::FGThruster(FGFDMExec *FDMExec, Element *el, int num ): FGForce(FDMEx PropertyManager->Tie( property_name, (FGForce *)this, &FGForce::GetYaw, &FGForce::SetYaw); if (el->GetName() == "direct") // this is a direct thruster. At this time - // only a direct thruster can be reversed. + // only a direct thruster can be reversed. { snprintf(property_name, 80, "propulsion/engine[%d]/reverser-angle-rad", EngineNum); PropertyManager->Tie( property_name, (FGThruster *)this, &FGThruster::GetReverserAngle, diff --git a/src/FDM/JSBSim/models/propulsion/FGThruster.h b/src/FDM/JSBSim/models/propulsion/FGThruster.h index c2adf1998..d37a73ad8 100644 --- a/src/FDM/JSBSim/models/propulsion/FGThruster.h +++ b/src/FDM/JSBSim/models/propulsion/FGThruster.h @@ -99,7 +99,7 @@ public: virtual void SetRPM(double rpm) {}; virtual double GetPowerRequired(void) {return 0.0;} virtual void SetdeltaT(double dt) {deltaT = dt;} - double GetThrust(void) {return Thrust;} + double GetThrust(void) const {return Thrust;} eType GetType(void) {return Type;} string GetName(void) {return Name;} void SetReverserAngle(double angle) {ReverserAngle = angle;} diff --git a/src/FDM/JSBSim/models/propulsion/FGTurbine.cpp b/src/FDM/JSBSim/models/propulsion/FGTurbine.cpp index 022d39235..16bc0fe1e 100644 --- a/src/FDM/JSBSim/models/propulsion/FGTurbine.cpp +++ b/src/FDM/JSBSim/models/propulsion/FGTurbine.cpp @@ -107,6 +107,8 @@ void FGTurbine::ResetToIC(void) double FGTurbine::Calculate(void) { + double thrust; + TAT = (Auxiliary->GetTotalTemperature() - 491.69) * 0.5555556; dt = State->Getdt() * Propulsion->GetRate(); ThrottlePos = FCS->GetThrottlePos(EngineNumber); @@ -144,19 +146,19 @@ double FGTurbine::Calculate(void) if (Seized) phase = tpSeize; switch (phase) { - case tpOff: Thrust = Off(); break; - case tpRun: Thrust = Run(); break; - case tpSpinUp: Thrust = SpinUp(); break; - case tpStart: Thrust = Start(); break; - case tpStall: Thrust = Stall(); break; - case tpSeize: Thrust = Seize(); break; - case tpTrim: Thrust = Trim(); break; - default: Thrust = Off(); + case tpOff: thrust = Off(); break; + case tpRun: thrust = Run(); break; + case tpSpinUp: thrust = SpinUp(); break; + case tpStart: thrust = Start(); break; + case tpStall: thrust = Stall(); break; + case tpSeize: thrust = Seize(); break; + case tpTrim: thrust = Trim(); break; + default: thrust = Off(); } - Thrust = Thruster->Calculate(Thrust); // allow thruster to modify thrust (i.e. reversing) + thrust = Thruster->Calculate(thrust); // allow thruster to modify thrust (i.e. reversing) - return Thrust; + return thrust; } //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% diff --git a/src/FDM/JSBSim/models/propulsion/FGTurbine.h b/src/FDM/JSBSim/models/propulsion/FGTurbine.h index 52193db99..c0524ce1b 100644 --- a/src/FDM/JSBSim/models/propulsion/FGTurbine.h +++ b/src/FDM/JSBSim/models/propulsion/FGTurbine.h @@ -169,7 +169,6 @@ public: double Calculate(void); double CalcFuelNeed(void); double GetPowerAvailable(void); -// double GetThrust(void) const {return Thrust;} /** A lag filter. Used to control the rate at which values are allowed to change. @param var a pointer to a variable of type double