Balloons behaving oddly at high altitudes

[jAllyn-ACTA]

I have been working with a high-altitude balloon model but it has been behaving very oddly when it is at high altitudes. All of the balloon models I have used work well at low altitudes, but at altitudes above ~60k feet, the balloons stop reacting correctly to winds. I have tried using my own model as well as the provided weather balloon model:
https://github.com/JSBSim-Team/jsbsim/blob/master/aircraft/weather-balloon/weather-balloon.xml

The balloons move at a constant velocity when the winds are constant (as expected), but they do not actually follow the winds. The balloons will travel only straight north or straight east. I am inputting the winds using the following code block:

      <run start="0" end="1" dt="0.1">

        <event name="Winds" continuous="true">

          <condition>simulation/sim-time-sec ge 0.0</condition>

          <set name="atmosphere/wind-north-fps" value="4.7983" />

          <set name="atmosphere/wind-east-fps" value="32.9567" />

          <set name="atmosphere/wind-down-fps" value="0" />

       </event>

      </run>

I am unsure if this is a bug with JSBSim or if I am doing something wrong. Any help would be greatly appreciated!

Some extra information:
I am running Octave to set up the winds depending on altitude. My octave code sets up the initial conditions, calls JSBSim, captures the output, and sets up the initial conditions for a new JSBSim run. The winds are calculated using GFS forecast data. Balloon path information is saved into a CSV file.

[dpculp]

I suspect this may be coming from the turbulence code.  I see in FGWinds.cpp that the turbulence is set to ttMilspec when I think it should be ttNone. That's my first guess.  Anyone else have another guess? Dave

…

On 6/3/22 11:32, jAllyn-ACTA wrote: I have been working with a high-altitude balloon model but it has been behaving very oddly when it is at high altitudes. All of the balloon models I have used work well at low altitudes, but at altitudes above ~60k feet, the balloons stop reacting correctly to winds. I have tried using my own model as well as the provided weather balloon model: https://github.com/JSBSim-Team/jsbsim/blob/master/aircraft/weather-balloon/weather-balloon.xml The balloons move at a constant velocity when the winds are constant (as expected), but they do not actually follow the winds. The balloons will travel only straight north or straight east. I am inputting the winds using the following code block: ` |<run start="0" end="1" dt="0.1"> <event name="Winds" continuous="true"> <condition>simulation/sim-time-sec ge 0.0</condition> <set name="atmosphere/wind-north-fps" value="4.7983" /> <set name="atmosphere/wind-east-fps" value="32.9567" /> <set name="atmosphere/wind-down-fps" value="0" /> </event> </run> | ` I am unsure if this is a bug with JSBSim or if I am doing something wrong. Any help would be greatly appreciated! Some extra information: I am running Octave to set up the winds depending on altitude. My octave code sets up the initial conditions, calls JSBSim, captures the output, and sets up the initial conditions for a new JSBSim run. The winds are calculated using GFS forecast data. Balloon path information is saved into a CSV file. — Reply to this email directly, view it on GitHub <#640>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/AFDBF2OTVKTVZYURXPDZSYTVNJFUBANCNFSM5XZZZGRA>. You are receiving this because you are subscribed to this thread.Message ID: ***@***.***>

[ermarch]

This, at least, explains why I get turbulence even on a wind free day in FlightGear.

[bcoconni]

Well, probability_of_exceedence_index is initialized to 0:

jsbsim/src/models/atmosphere/FGWinds.cpp

Lines 91 to 93 in 59389f6

	// Milspec turbulence model
	windspeed_at_20ft = 0.;
	probability_of_exceedence_index = 0;

and ttMilSpec returns no turbulence wind when probability_of_exceedence_index is equal to 0:

jsbsim/src/models/atmosphere/FGWinds.cpp

Lines 246 to 255 in 59389f6

	case ttMilspec:
	case ttTustin: {
	// an index of zero means turbulence is disabled
	// airspeed occurs as divisor in the code below
	if (probability_of_exceedence_index == 0 || in.V == 0) {
	vTurbulenceNED(eNorth) = vTurbulenceNED(eEast) = vTurbulenceNED(eDown) = 0.0;
	vTurbPQR(eP) = vTurbPQR(eQ) = vTurbPQR(eR) = 0.0;
	return;
	}

So the ttMilSpec turbulence code is sort of disabled by default, unless one set the property atmosphere/turbulence/milspec/severity to a value different than 0.

jsbsim/src/models/atmosphere/FGWinds.cpp

Lines 557 to 559 in 59389f6

	PropertyManager->Tie("atmosphere/turbulence/milspec/severity",
	this, &FGWinds::GetProbabilityOfExceedence,
	&FGWinds::SetProbabilityOfExceedence);

jsbsim/src/models/atmosphere/FGWinds.h

Lines 269 to 271 in 59389f6

	/// allowable range: 0-7, 3=light, 4=moderate, 6=severe turbulence
	virtual void SetProbabilityOfExceedence( int idx) {probability_of_exceedence_index = idx;}
	virtual int GetProbabilityOfExceedence() const { return probability_of_exceedence_index;}

I agree that setting turbulence to ttNone would make things clearer.

[bcoconni]

Anyway @jAllyn-ACTA, it's easy to check if there is some turbulence wind, just check if properties atmosphere/turb-north-fps, "atmosphere/turb-east-fps and "atmosphere/turb-down-fps are non-zero.

[dpculp]

I just made a test flight through 75K feet with an aircraft model and the wind direction and speed did not vary from the initial conditions (70 degrees, 40 knots for this test). It looks like the wind changes are not coming from the atmosphere or wind models.

I will next try a balloon model, but I have no experience doing this, so it may take some time to figure it out. It could be that the problem is in the balloon code - how it reacts to the wind might have a discontinuity at about 60K feet?

@bcoconni Sorry about the noise. I think initializing with the ttMilspec turbulence model is fine as long as there is no way to "accidentally" have it supply perturbations, and you're right - there doesn't seem to be a way to accidentally do this.

@ermarch In Flightgear you're probably using the the Flightgear atmosphere and weather system (the default) and last I looked there are some non-standard-atmosphere effects that you can't turn off.

[bcoconni]

@bcoconni Sorry about the noise. I think initializing with the ttMilspec turbulence model is fine as long as there is no way to "accidentally" have it supply perturbations, and you're right - there doesn't seem to be a way to accidentally do this.

No problem. This was not obvious to me either, but I was surprised that such a bug would not have been reported for years. Following my findings, I was relieved that the code was sound 😄

[dpculp]

I ran a test using the "weather-balloon" model that is included with the JSBSim aircraft directory. I edited the initialization file, reset10.xml, to include some wind:

winddir 70.0
vwind 40.0

(sorry, can't get the actual code to appear correctly here)

I changed the script, scripts/weather-balloon.xml, to report the ground track and ground speed:

    <notify>
        <property>position/h-agl-ft</property> 
        <property>velocities/vg-fps</property>
        <property>flight-path/psi-gt-rad</property>
        <property>metrics/radius-ft</property>
    </notify>

I also changed the script to terminate when exceeding 80K feet altitude.

The ground speed and track held steady all the way. I'm not able to duplicate the error you're seeing.

[bcoconni]

@dpculp, I edited your comment with the triple back-tick formatting as @seanmcleod mentioned below.

[dpculp]

Thanks, Bertrand.

[jAllyn-ACTA]

@dpculp I was seeing the issue only when I started the balloon from high altitude, not when I start it from the ground. If I start from the ground, everything works fine. If I start at higher altitudes, the issue occurs.

[seanmcleod]

(sorry, can't get the actual code to appear correctly here)

@dpculp take a look at the triple back-tick formatting options and the syntax highlighting option - https://docs.github.com/en/get-started/writing-on-github/working-with-advanced-formatting/creating-and-highlighting-code-blocks

I typically use it with python, c++ and sh (for showing output from say JSBSim executable, showing a command line example etc.).

[bcoconni]

I was also wondering if the air density would not be too thin above 60kft for the wind to produce a significant drag ? I mean the weather balloon example uses the good old $\bar{q}=\frac{1}{2}\rho V^2$ aka aero/qbar-psf:

jsbsim/aircraft/weather-balloon/weather-balloon.xml

Lines 217 to 238 in 59389f6

	<function name="aero/force/BasicDrag">
	<description>Drag. Table from [1]</description>
	<product>
	<value name="Tweak"> 1.22 </value>
	<property>aero/qbar-psf</property>
	<property>metrics/Sref-sqft</property>
	<table>
	<independentVar>aero/Re-cor</independentVar>
	<tableData>
	8766 0.39
	28687 0.44
	81039 0.50
	110575 0.49
	232685 0.42
	460735 0.31
	613333 0.24
	690250 0.21
	883256 0.20
	</tableData>
	</table>
	</product>
	</function>

[seanmcleod]

And in addition to triple back-tick code formatting we now also have Github's support for Latex math! 😉

[jAllyn-ACTA]

@bcoconni Even if there isn't significant drag, the balloon's horizontal velocity should not exceed the wind's horizontal velocity. If I start the balloon at a high altitude, the problem occurs. Starting from the ground and going up to these altitudes does not cause this issue but starting from a high altitude does. I need to start at high altitudes for the work I am trying to do.

Another issue I found was that the balloon velocity seems to be reasonable (although in the wrong direction) until either the NS or EW wind magnitude shifts to negative. Am I using the wind block in my original post correctly? I am still new to JSBSim so I could be doing something wrong.

Thanks for your help!

[seanmcleod]

@jAllyn-ACTA could you post the contents of the initialization file that you're using with the default balloon model in JSBSim which demonstrates the issue you're seeing, with the relevant properties being logged which show the unexpected velocities you're seeing versus what you would expect.

Which would then make it easier for us to reproduce the issue, otherwise a fair bit of time is wasted trying to setup a test to try and reproduce what you're seeing and often fail to capture the exact requirements that are needed to reproduce the issue.

[jAllyn-ACTA]

Initialization file:

<?xml version="1.0" encoding="utf-8"?>
<initialize name="cruise">
	<altitude unit="FT">100000.0</altitude>
	<vt unit="FT/SEC">0</vt>
	<gamma unit="DEG">0</gamma>
	<latitude unit="DEG">36.9297</latitude>
	<longitude unit="DEG">-111.4489</longitude>
	<phi unit="DEG">0</phi>
	<psi unit="DEG">0</psi>
	<beta unit="DEG">0</beta>
</initialize>

Script file:

<?xml version="1.0" encoding="utf-8"?>
<?xml-stylesheet type="text/xsl" href="http://jsbsim.sf.net/JSBSimScript.xsl"?>
<runscript xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance"
	xsi:noNamespaceSchemaLocation="http://jsbsim.sf.net/JSBSimScript.xsd"
	name="Balloon Flight">
<use aircraft="weather-balloon" initialize="balloonInit"/>
<run start="0" end="1" dt="0.01">
	<event name="Winds" continuous="true">
		<condition>simulation/sim-time-sec ge 0.0</condition>
		<set name="atmosphere/wind-north-fps" value="0" />
		<set name="atmosphere/wind-east-fps" value="0" />
		<set name="atmosphere/wind-down-fps" value="0" />
	</event>
</run>
	<output name="Output/example_output.csv" type="CSV" rate="100" file="../scripts/unitconversions.xml">
		<property caption="Latitude"> position/lat-gc-deg </property>
		<property caption="Longitude"> position/long-gc-deg </property>
		<property caption="Altitude (sea level)">  position/h-sl-ft </property> <!-- height above sea level -->
		<property caption="Altitude (ground level)"> position/h-agl-ft </property> <!-- height above ground level -->
		<property caption="Airspeed (ft/s)">velocities/vt-fps</property>
		<property apply="convert-rad-To-deg" caption="Phi (deg)">attitude/phi-rad</property>
		<property apply="convert-rad-To-deg" caption="Theta (deg)">attitude/theta-rad</property>
		<property apply="convert-rad-To-deg" caption="Psi (deg)">attitude/psi-rad</property>
		<property caption="Beta (deg)">ic/beta-deg</property>
		<property caption="Gamma (deg)">ic/gamma-deg</property>
		<property caption="Qbar (psf)">aero/qbar-psf</property>
	</output>
</runscript>

Example output file:
example_output.csv

[jAllyn-ACTA]

@bcoconni I have checked multiple times to make sure that all of the values transfer between Octave and JSBSim correctly. The Initialization and Script files are modified by Octave every "Octave timestep." The "Octave timesteps" include multiple "JSBSim timesteps." Here is the "main loop" where Octave calls JSBSim:

% Main Loop
  for i=1:totalDt
    
   % Find altitude bar
    index = 1;
    for j=1:length(windHeight)
      if (windHeight(j) > altitude)
        break;
      endif
      index = j;
    endfor
    
    % Save current winds
    meanWindNS = currWindNS(index);
    meanWindEW = currWindEW(index);
    
    % Add realizations to winds
    usedWindNS = meanWindNS + realizationNS(index);
    usedWindEW = meanWindEW + realizationEW(index);
    
    % Write Initial Conditions File
    fid = fopen(tempInit_xml, 'w');
    fprintf(fid, CreateInitCondXMLStr(altitude, velocity, gamma, latitude, longitude, phi, psi, beta));
    fclose(fid);
    
    % Write Script File
    fid = fopen(tempScript_xml, 'w');
    fprintf(fid, CreateScriptXMLStr(output_csv, dt, ndt, usedWindNS, usedWindEW, currWindV, aircraftName));
    fclose(fid);

    % Run JSBSim
    system(['C:/Jake/Software/JSBSim/JSBSim.exe ' tempScript_xml]);
    
    % Read CSV
    outCSV = dlmread(output_csv,',',1,0);
    
    % Set new ICs
    row = length(outCSV(:,1));
    latitude = outCSV(row, 2);
    longitude = outCSV(row, 3);
    altitude = outCSV(row, 4);
    velocity = outCSV(row, 6);
    phi = outCSV(row, 7);
    psi = outCSV(row, 9);
    beta = outCSV(row, 10);
    gamma = outCSV(row, 11);
    qbar = outCSV(row, 12);
    
    % Write to CSV
    fprintf(writeCSV, "%d, %d, %d, %d, %s, %d, %d, %d, %d\n",... 
    mc, latitude, longitude, altitude, datestr(currTime, "dd-mmm-yyyy HH:MM:SS.FFF"),... 
    usedWindNS, usedWindEW,... %meanWindNS, meanWindEW,... %sigmaNS, sigmaEW,... 
    velocity, qbar);%, phi, psi, beta, gamma);
    
    % Increment time
    currTime = datevec(addtodate(datenum(currTime), timeDeltaPerCall, "second"));
    
    % End simulation on ground impact
    if(endOnGroundImpact == "yes" && altitude <= 0.1)
      break;
    endif
  endfor

Another option I thought of is to create a table of winds based on altitude and run JSBSim all at once instead of iterating this way. Is it possible to set up the script file to have a table based on altitude (like the balloon's drag table)? Doing it this way may help determine if my Octave script is causing the issue or not.

[bcoconni]

Are you executing JSBSim.exe just as many time as time steps you want to simulate ? i.e. you run the simulation for a single time step $\delta t$ ?

If that is so, then this is no surprise that you get silly results at least because multistep iteration algorithms (Adams-Bashforth) are not designed for that. JSBSim always initializes the multistep algorithm assuming that the vehicle was in a steady state for $t&lt;0$ so this basically corrupts the balloon dynamics if you run a fresh simulation at each time step.

Another option I thought of is to create a table of winds based on altitude and run JSBSim all at once instead of iterating this way. Is it possible to set up the script file to have a table based on altitude (like the balloon's drag table)?

Yes, check the script scripts/J2461.xml for an example:

jsbsim/scripts/J2461.xml

Lines 21 to 78 in 4ab60d2

	<event name="Winds aloft" continuous="true">
	<condition>simulation/sim-time-sec ge 10000.05</condition>
	<set name="atmosphere/psiw-rad">
	<function>
	<table>
	<independentVar lookup="row"> position/h-agl-ft </independentVar>
	<tableData>
	<!-- Alt (ft) Dir (rad) -->
	0.0 0.80
	5000.0 0.90
	10000.0 1.00
	15000.0 0.90
	20000.0 0.70
	25000.0 0.80
	30000.0 0.90
	35000.0 1.00
	40000.0 1.10
	45000.0 1.20
	50000.0 1.30
	55000.0 1.40
	60000.0 1.30
	65000.0 1.20
	70000.0 1.30
	100000.0 1.40
	200000.0 1.50
	</tableData>
	</table>
	</function>
	</set>
	<set name="atmosphere/wind-mag-fps">
	<function>
	<table>
	<independentVar lookup="row"> position/h-agl-ft </independentVar>
	<tableData>
	<!-- Alt (ft) Speed (fps) -->
	0.0 10.0
	5000.0 40.0
	10000.0 60.0
	15000.0 120.0
	20000.0 135.0
	25000.0 140.0
	30000.0 130.0
	35000.0 120.0
	40000.0 110.0
	45000.0 110.0
	50000.0 100.0
	55000.0 90.0
	60000.0 85.0
	65000.0 60.0
	70000.0 40.0
	100000.0 10.0
	200000.0 5.0
	</tableData>
	</table>
	</function>
	</set>
	<notify/>
	</event>

[bcoconni]

If that is so, then this is no surprise that you get silly results at least because the multistep iteration algorithms (Adams-Bashforth) are not designed for that. JSBSim always initializes the multistep algorithm assuming that the vehicle was in a steady state for so this basically corrupts the balloon dynamics if you run a fresh simulation at each time step.

Also JSBSim does not provide a means to initialize the accelerations so your initial accelerations are most likely (very) different than what they were at the previous time step.

[jAllyn-ACTA]

OK, thanks for the response. I will try the table method and see if it fixes the issue.

Are you executing JSBSim.exe just as many time as time steps you want to simulate ? i.e. you run the simulation for a single time step dt?

Not quite, it runs about 100 timesteps in JSBSim before resetting in Octave. Hopefully doing the table lookup for the winds will fix this!

[bcoconni]

OK, thanks for the response. I will try the table method and see if it fixes the issue.

Well, if it doesn't at least we will have a reproducible test case to check what is going wrong 😉

[seanmcleod]

Hmm, does seem like a fairly convoluted setup now that you've shared more details. I see in the meantime @bcoconni has made some suggestions in terms of wind table lookups for a single JSBSim run. When I initially saw your Octave JSBSim setup with multiple runs of JSBSim to handle updating the wind at different altitudes as the balloon's altitude changed my thinking was an option along these lines.

import jsbsim

fdm = jsbsim.FGFDMExec() 
fdm.load_model('weather-balloon')

fdm['ic/h-sl-ft'] = 100000
fdm.run_ic()

while fdm['position/h-sl-ft'] > 0.1:
  windN, windE, windD = lookupWind(fdm['position/h-sl-ft'])
  fdm['atmosphere/wind-north-fps'] = windN
  fdm['atmosphere/wind-east-fps'] = windE
  fdm['atmosphere/wind-down-fps'] = windD

   fdm.run()

   # Log state to csv file

[jAllyn-ACTA]

@bcoconni @seanmcleod Thank you both for your help and suggestions. Using the wind table method solved the issue. The weather balloon model is performing as expected now.