JSBSim-Team/jsbsim

FGStandardAtmosphere:CalculatePressureAltitude() can cause NaN values

vranki opened this issue · 19 comments

I'm submitting a ...

  • bug report
  • feature request
  • support request => Please do not submit support request here, see note at the top of this template.

Describe the issue
When studying #809 I tracked a source for NaN values causing jsbsim to crash.
FGStandardAtmosphere:CalculatePressureAltitude() with pressure = 0 causes pressure_altitude to be set to log(0). This results in NaNs in future calculations.

https://github.com/JSBSim-Team/jsbsim/blob/master/src/models/atmosphere/FGStandardAtmosphere.cpp#L549

What is the current behavior?
Situation when pressure = 0 is not handled correctly and causes jsbsim to crash.

What is the expected behavior?
Return valid value (or error/assertion, if not possible).

What is the motivation / use case for changing the behavior?
Not crashing.

Please tell us about your environment:

  • OS Ubuntu 20.04
  • JSBSim version 1.1.13
  • Language if applicable C++

Other information

Getting into situation where pressure = 0 is probably caused by another bug elsewhere, but this should be handled correctly anyway.

I just checked to see if flying below ground level could cause a bad pressure value, and I got good values as far as 12,000 feet below ground level, which is as far as I went.

FYI

Try going to space, maybe that'll trigger the issue. I believe I got it after warping to something like 3000000ft agl due to unrelated bug.

There is a Space Shuttle model for JSBSim/FlightGear that has been flying to space for years without issues 😉

I don't have time for this.

Interesting but is there a realistic environment/chance of reaching a 0.0 pressure? Outer space is 0.00000000001322 Pa
Also for my own curiosity I tested JSBSim limits of altitude above earth. And I reached 38,834,063 ft above seal level before the simulation crashed. That's 7,354 miles.
For reference the International Space Station is 260 miles above earth and the Space Shuttles' maximum operating orbit was 390 miles, but never flew beyond 360 miles from earth.

I'm not sure handling of a 0.0 pressure, even with a check/asset, would solve anything except in fringe cases caused by other factors, which should probably be addressed separately?

Ok. Perhaps zero value was caused by some earlier miscalculation somewhere. I'd suggest adding assert(pressure > 0) to stop simulation on invalid input value.

Try going to space, maybe that'll trigger the issue. I believe I got it after warping to something like 3000000ft agl due to unrelated bug.

This kind of mishap can happen when the vehicle hits violently the ground and compresses its gears and contact points to a ridiculous amount. This can sometimes result in propelling the vehicle to outer space which is of course completely unrealistic.

Ok. Perhaps zero value was caused by some earlier miscalculation somewhere. I'd suggest adding assert(pressure > 0) to stop simulation on invalid input value.

To my humble point of view, asserts are meant to never being triggered. Each time I have added an assert in JSBSim it was to check that an event would never occur. And would that event occur this would mean that there is an error in the code and that error should be fixed.

Regarding the current issue, we definitely know that the event P==0.0 can occur so we should either make sure that this event should never happen or we couldreport gracefully the error to the user. I'm pretty sure asserts do not fall in either category.

@gallonmate remark about pressure never dropping to the absolute zero is interesting. Indeed JSBSim atmosphere is only truthful to the ISA 1976 model up to 91 km (approx. 300,000ft), beyond that extrapolation takes place and I am not surprised that at some point FGStandardAtmosphere returns garbage. JSBSim is not (yet ?) designed for interplanetary travel 😉.

// GeoPot Alt Temp GeoPot Alt GeoMet Alt
// (ft) (deg R) (km) (km)
// ----------- -------- ---------- ----------
StdAtmosTemperatureTable << 0.0000 << 518.67 // 0.000 0.000
<< 36089.2388 << 389.97 // 11.000 11.019
<< 65616.7979 << 389.97 // 20.000 20.063
<< 104986.8766 << 411.57 // 32.000 32.162
<< 154199.4751 << 487.17 // 47.000 47.350
<< 167322.8346 << 487.17 // 51.000 51.413
<< 232939.6325 << 386.37 // 71.000 71.802
<< 278385.8268 << 336.5028 // 84.852 86.000
<< 298556.4304 << 336.5028; // 91.000 - First layer in high altitude regime

All this needs some further investigation in order to address the issue that you have reported @vranki .

JSBSim is not (yet ?) designed for interplanetary travel

Yep, we have a Mars atmosphere available, but as you mentioned nothing in-between 91km above earth, so the vehicle would have to teleport 😉

In the discussion #809 (comment) @vranki mentioned that due to some time roll over bug he was setting delta-T to some invalid temperature offset of -4000 which he thought was then 'moving the aircraft to outer space'.

What would be the correct return value if pressure was zero? We could just return that in the situation and avoid log(0). That would make interplanetary travel possible in future.

What would be the correct return value if pressure was zero?

Well, assuming there is a point where the pressure really goes to 0, then there are infinite altitudes greater than that which are valid answers for a pressure of 0.

In the discussion #809 (comment) @vranki mentioned that due to some time roll over bug he was setting delta-T to some invalid temperature offset of -4000 which he thought was then 'moving the aircraft to outer space'.

Setting delta-T to -4000 is most certainly a path to disaster since that would set the ambient temperature to values below the absolute zero (i.e. negative temperatures on the Rankine scale) which is a physical nonsense. For example that would result in pressure increasing as the vehicle altitude is increasing.

This is certainly an error that should be handled by JSBSim (currently it is not).

What would be the correct return value if pressure was zero?

Well, assuming there is a point where the pressure really goes to 0, then there are infinite altitudes greater than that which are valid answers for a pressure of 0.

Given that LapseRates[7]==0.0, the pressure at altitudes above 91km are computed using the exponential formula just below:

return PressureBreakpoints[b]*exp(-g0*deltaH/(Rdry*Tmb));

According to a quick test I ran on Python (under Windows 11 - 64 bits), exp(-800.0)==0.0. In other words, the exponential function decreases so fast that exp(-800.0) cannot be correctly represented in double precision and the FPU rounds its value down to 0.0. I presume that we'd get a similar result in C++ (I have not tested).

This means that as soon as g0*deltaH/(Rdry*Tmb)>=800.0 the pressure is set to zero due to a floating point underflow and this could occur either when deltaH is sufficiently high (such as in @gallonmate's test: #815 (comment)) or when the ambient temperature Tmb is getting sufficiently close to 0.0.

I think this boils down to that in my opinion JSBSim should try to detect invalid values and do whatever is needed to return false from FDMExec->Run(); to halt the simulation. Invalid values can also be generated in hard crash, so it's not (only) API misuse.

Current functionality of just ignoring invalid values until some assert causes the whole process to crash is not good IMO. Returning false from Run() would be a clean solution when the simulation cannot continue.

Some cases are easier than others. For example a delta-T offset which then results in temperatures below absolute zero is an easier case to check and handle.

Pressure reaching zero at very high altitudes is a bit trickier. Since the ISA atmosphere model tops out around 91km, do we simply clip past this point, or do we throw an exception, or allow interpolation to continue to ~7,354 miles as per @gallonmate's calculation even though the ISA atmosphere model doesn't cover this region, and then throw an exception at that point.

Throwing an exception at 91km would break things like the Space Shuttle FlightGear model.

In terms of inputs being out of range, remember that for FDM data in terms of lookup tables, e.g. coefficient of lift versus alpha, that JSBSim clips the output, it doesn't extrapolate and doesn't throw an exception.

Talking of the atmosphere above 91km specifically, the other option is to implement one of these models above 91km. Although they also have limits, e.g. 2000km for JB2008 so we would still need to decide what to do when passing their limits.

https://en.wikipedia.org/wiki/International_Standard_Atmosphere

NRLMSISE-00 is a newer model of the Earth's atmosphere from ground to space, developed by the US Naval Research Laboratory taking actual satellite drag data into account. A primary use of this model is to aid predictions of satellite orbital decay due to atmospheric drag. The COSPAR International Reference Atmosphere (CIRA) 2012 and the ISO 14222 Earth Atmosphere Density standard both recommend NRLMSISE-00 for composition uses.

JB2008 is a newer model of the Earth’s atmosphere from 120 km to 2000 km, developed by the US Air Force Space Command and Space Environment Technologies taking into account realistic solar irradiances and time evolution of geomagnetic storms.[13] It is most useful for calculating satellite orbital decay due to atmospheric drag. Both CIRA 2012 and ISO 14222 recommend JB2008 for mass density in drag uses.

Current functionality of just ignoring invalid values until some assert causes the whole process to crash is not good IMO. Returning false from Run() would be a clean solution when the simulation cannot continue.

@vranki just to avoid a few misunderstandings:

  • The fact that JSBSim does not handle some errors (including the errors that you have reported in this issue) is not made on purpose. We are trying to address errors as they are discovered, knowing that not all of them have a simple way of being detected, not to mention about finding a way to fix or avoid them.
  • I was mentioning previously that asserts are not a good way to address the current issue. However by no means was I meaning to keep the current situation unmodified. If the plan was to ignore these errors, this issue would have already been closed.

Some cases are easier than others. For example a delta-T offset which then results in temperatures below absolute zero is an easier case to check and handle.

Regarding this issue, I guess that the simpler way to fix it is to make sure that the temperature is always above zero. Something along these lines:

  T += TemperatureBias;

  if (GeoPotAlt <= GradientFadeoutAltitude)
    T += TemperatureDeltaGradient * GradientFadeoutAltitude;

-  return T;
+  return max(T, 1.0);

However how should we report this:

  • Should JSBSim swallow the error and continue execution as if nothing happened ?
  • Should a single error message be reported in the log messages and execution continue ?
  • Should JSBSim throw an exception ? This option will most likely trigger the termination of the calling program ?

The latter option is certainly the safer from our point of view but it might be frustrating for the user as it might not be obvious to relate a below-absolute-zero-temperature-error to the source of the problem i.e. aggressive values being set in atmosphere/delta-T or atmosphere/SL-graded-delta-T.

Pressure reaching zero at very high altitudes is a bit trickier. Since the ISA atmosphere model tops out around 91km,

The ISA atmosphere is not limited to 91km, it provides a model up to 1000 km (see for example the introduction to §1.3)

The thing is that JSBSim grossly extrapolates atmosphere properties above 86km because the ISA model is getting significantly more complex above that altitude. Above that threshold, species should be computed separately instead of addressing a single homogeneous gas (the air).

Nothing infeasible but it is definitely more involved than the current pressure exponential damping 😄

Talking of the atmosphere above 91km specifically, the other option is to implement one of these models above 91km. Although they also have limits, e.g. 2000km for JB2008 so we would still need to decide what to do when passing their limits.

Not sure how complex these models are but would they be used, I guess we should consider using a library for such a level of sophistication. Just as GeographicLib has been introduced because at some point the code for geodetic computations is getting so specialized that it deserves a library of its own.

I just noticed there already is NRLMSISE-00 model in the atmosphere folder. JB2008 is very similar and more accurate but it's table data is only for 120km to 2000km. The two models could be merged or the basic data could be extracted and used to extend the tables in Standard Atmosphere.cpp

/** Models the MSIS-00 atmosphere.
This is a wrapper for the NRL-MSIS-00 model 2001:
This C++ format model wraps the NRLMSISE-00 C source code package - release
20020503

But as @seanmcleod mentioned, none of these models would solve the NaN pressure issue. Both models use extrapolation in their final table entry and will likely cause NaN values at some very far distances from earth.

My suggestion is simply put a clamp in StandardAtmosphere.cpp's final extrapolation, so that pressure can never reach NaN or zero. Outer Space pressure is regarded as 1.322 × 10-11 Pa. I'd use this number. Maybe send a message to Log if this number is reached. "Vehicle has reached altitude with almost zero pressure found in Outer Space. If this was unintentional please check all inputs for valid numbers or this vehicle possibly crashed hard into the ground, launching it into Outer Space". 😃

Although this does not solve the issue of the user manually initializing a ground level pressure of zero, which also causes JSBSim to crash.

Here is a general list of pressures, with the bottom number showing 500km above Earth, all the way up to lowest theoretical pressure. I guess a decision would need to be made if zero pressure should ever be able to simulate in JSBSim or if clamping the entire pressure in JSBSim to the lowest non-zero pressure that it's double precision can handle.
https://en.wikipedia.org/wiki/Orders_of_magnitude_(pressure)
image

@bcoconni Brings up another good point, Temperature should not and cannot go below absolute zero. This is actually noted and handled in the Mars Atmosphere.

// LIMIT the temperatures so they do not descend below absolute zero.

In summary:

  • zero and NaN pressure needs to be handled. Issues caused by StandardAtmosphere.cpp's calculation and also by user setting or initializing zero and low pressures.
  • temperatures need to be handled to not go below absolute zero.

Suggested solutions to discuss:

  • Clamp all functions related to Set Pressure, so that they can never be zero. Clamp StandardAtmosphere.cpp's final pressure so it never reaches zero/NaN. Send a message to Log.
  • Clamp all Set Temperature functions so that Temperature cannot go below absolute zero. Send a message to log.

Suggested solutions to discuss:

  • Clamp all functions related to Set Pressure, so that they can never be zero. Clamp StandardAtmosphere.cpp's final pressure so it never reaches zero/NaN. Send a message to Log.
  • Clamp all Set Temperature functions so that Temperature cannot go below absolute zero. Send a message to log.

This has been implemented by the PR #836. This issue is now considered fixed.