1 /*****************************************************************************
4 Author: Christian Mayer
8 -------- Copyright (C) 1999 Christian Mayer (fgfs@christianmayer.de) --------
10 This program is free software; you can redistribute it and/or modify it under
11 the terms of the GNU General Public License as published by the Free Software
12 Foundation; either version 2 of the License, or (at your option) any later
15 This program is distributed in the hope that it will be useful, but WITHOUT
16 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
17 FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
20 You should have received a copy of the GNU General Public License
21 along with this program; if not, write to the Free Software
22 Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
24 Further information about the GNU General Public License can also be found on
25 the world wide web at http://www.gnu.org.
27 FUNCTIONAL DESCRIPTION
28 ------------------------------------------------------------------------------
29 A hot air balloon simulator
32 ------------------------------------------------------------------------------
33 01.09.1999 Christian Mayer Created
34 03.10.1999 Christian Mayer cleaned the code by moveing WeatherDatabase
35 calls inside the update()
36 *****************************************************************************/
38 /****************************************************************************/
40 /****************************************************************************/
50 #include <simgear/constants.h>
52 #include <Aircraft/aircraft.hxx>
54 #include "BalloonSim.h"
58 /****************************************************************************/
59 /********************************** CODE ************************************/
60 /****************************************************************************/
62 /****************************************************************************/
65 /* Set the balloon model to some values that seem reasonable as I haven't */
66 /* got much original values */
68 /****************************************************************************/
72 ground_level = 3400.0;
74 sgSetVec3(gravity_vector, 0.0, 0.0, -9.81);
75 sgSetVec3(velocity, 0.0, 0.0, 0.0);
76 sgSetVec3(position, 0.0, 0.0, 0.0);
77 sgSetVec3(hpr, 0.0, 0.0, 0.0);
79 /************************************************************************/
80 /* My balloon has a radius of 8.8 metres as that gives a envelope */
81 /* volume of about 2854 m^3 which is about 77000 ft^3, a very common */
82 /* size for hot air balloons */
83 /************************************************************************/
85 balloon_envelope_area = 4.0 * (8.8 * 8.8) * SGD_PI;
86 balloon_envelope_volume = (4.0/3.0) * (8.8 * 8.8 * 8.8) * SGD_PI;
88 wind_facing_area_of_balloon = SGD_PI * (8.8 * 8.8);
89 wind_facing_area_of_basket = 2.0; //guessed: 2 m^2
91 cw_envelope=0.45; //a sphere in this case
94 weight_of_total_fuel = 40.0; //big guess
95 weight_of_envelope = 200.0; //big guess
96 weight_of_basket = 40.0; //big guess
97 weight_of_cargo = 750.0; //big guess
100 max_flow_of_fuel_per_second=10.0*1.0/3600.0; //assuming 10% of one hour of total burn time
101 current_burner_strength = 0.0; //the throttle
103 lambda = 0.15; //for plasic
104 l_of_the_envelope = 1.0/1000.0; //the thickness of the envelope (in m): 1mm
106 T = 273.16 + 130.6; //Temperature in the envelope => still at ground level
109 void balloon::update()
111 /************************************************************************/
112 /* I'm simplifying the balloon by reducing the simulation to two */
114 /* the center of the basket (CB) and the center of the envelope (CE) */
119 /* cg (=center of gravity) */
123 /* On each center are forces acting: gravity and wind resitance. CE */
124 /* additionally got the lift (=> I need to calculate the weight of the */
125 /* air inside, too) */
127 /* The weight of the air in the envelope is dependant of the tempera- */
128 /* ture. This temperature is decreasing over the time that is dependant */
129 /* of the insulation of the envelope material (lambda), the gas used */
130 /* (air) and the wind speed. For a plane surface it's for air: */
132 /* alpha = 4.8 + 3.4*v with v < 5.0 m/s */
134 /* The value k that takes all of that into account is defined as: */
136 /* 1 / k = 1 / alpha1 + 1 / alpha2 + l / lambda */
138 /* with 'l' beeing the 'length' i.e. thickness of the insulator, alpha1 */
139 /* the air inside and alpha2 the air outside of the envelope. So our k */
142 /* k = 1 / [1/4.8 + 1/(4.8+3.4v) + l/lambda] */
144 /* The energy lost by this process is defined as: */
146 /* dQ = k * A * t * dT */
148 /* with Q being the energy, k that value defined above, A the total */
149 /* area of the envelope, t the time (in hours) and dT the temperature */
150 /* difference between the inside and the outside. */
151 /* To get the temperature of the air in the inside I need the formula: */
153 /* dQ = cAir * m * dT */
155 /* with cAir being the specific heat capacity(?) of air (= 1.00 kcal / */
156 /* kg * degree), m the mass of the air and dT the temperature change. */
157 /* As the envelope is open I'm assuming that the same air pressure is */
158 /* inside and outside of it (practical there should be a slightly */
159 /* higher air pressure in the inside or the envelope would collapse). */
160 /* So it's easy to calculate the density of the air inside: */
162 /* rho = p / R * T */
164 /* with p being the pressure, R the gas constant(?) which is for air */
165 /* 287.14 N * m / kg * K and T the absolute temperature. */
167 /* The value returned by this function is the displacement of the CB */
168 /************************************************************************/
170 /************************************************************************/
171 /* NOTE: This is the simplified version: I'm assuming that the whole */
172 /* balloon consists only of the envelope.I will improove the simulation */
173 /* later, but currently was my main concern to get it going... */
174 /************************************************************************/
176 // I realy don't think there is a solution for this without WeatherCM
177 // but this is a hack, and it's working -- EMH
181 // gain of energy by heating:
182 if (fuel_left > 0.0) //but only with some fuel left ;-)
184 float fuel_burning = current_burner_strength * max_flow_of_fuel_per_second * dt * weight_of_total_fuel; //in kg
186 //convert to cubemetres (I'm wrongly assuming 'normal' conditions; but that's correct for my special case)
187 float cube_metres_burned = fuel_burning / 2.2; //2.2 is the density for propane
189 fuel_left -= fuel_burning / weight_of_total_fuel;
191 // get energy through burning.
192 Q += 22250.0 * cube_metres_burned; //22250 for propan, 29500 would be butane and if you dare: 2580 would be hydrogen...
195 // calculate the new temperature in the inside:
196 T += Q / (1.00 * mAir);
198 //calculate the masses of the envelope and the basket
199 float mEnvelope = mAir + weight_of_envelope;
200 float mBasket = weight_of_total_fuel*fuel_left + weight_of_basket + weight_of_cargo;
202 float mTotal = mEnvelope + mBasket;
204 //calulate the forces
205 sgVec3 fTotal, fFriction, fLift;
207 sgScaleVec3(fTotal, gravity_vector, mTotal);
209 //sgAddVec3(fTotal, fLift); //FIXME: uninitialized fLift
210 //sgAddVec3(fTotal, fFriction); //FIXME: uninitialized fFriction
212 //claculate acceleration: a = F / m
213 sgVec3 aTotal, vTotal, dTotal;
215 sgScaleVec3(aTotal, fTotal, 1.0 / mTotal);
217 //integrate the displacement: d = 0.5 * a * dt**2 + v * dt + d
218 sgScaleVec3(vTotal, velocity, dt);
219 sgScaleVec3(dTotal, aTotal, 0.5*dt*dt); sgAddVec3(dTotal, vTotal);
221 //integrate the velocity to 'velocity': v = a * dt + v
222 sgScaleVec3(vTotal, aTotal, dt); sgAddVec3(velocity, vTotal);
224 /************************************************************************/
225 /* VERY WRONG STUFF: it's just here to get some results to start with */
226 /************************************************************************/
228 // care for the ground
229 if (position[2] < (ground_level+0.001) )
230 position[2] = ground_level;
233 sgAddVec3(position, dTotal);
235 //cout << "BallonSim: T: " << (T-273.16) << " alt: " << position[2] << " ground: " << ground_level << " throttle: " << current_burner_strength << "\n";
238 void balloon::set_burner_strength(const float bs)
240 if ((bs>=0.0) && (bs<=1.0))
241 current_burner_strength = bs;
244 void balloon::getVelocity(sgVec3 v) const
246 sgCopyVec3(v, velocity);
249 void balloon::setVelocity(const sgVec3 v)
251 sgCopyVec3(velocity, v);
254 void balloon::getPosition(sgVec3 v) const
256 sgCopyVec3(v, position);
259 void balloon::setPosition(const sgVec3 v)
261 sgCopyVec3(position, v);
264 void balloon::getHPR(sgVec3 angles) const //the balloon isn't allways exactly vertical
266 sgCopyVec3(angles, hpr);
269 void balloon::setHPR(const sgVec3 angles) //the balloon isn't allways exactly vertical
271 sgCopyVec3(hpr, angles);
274 void balloon::setGroundLevel(const float altitude)
276 ground_level = altitude;
279 float balloon::getTemperature(void) const
284 float balloon::getFuelLeft(void) const
295 sgVec3 pos={0.0, 0.0, 0.0};
300 bool hysteresis = false; // moving up
303 for (int i=0; i<100; i++)
305 bal.update(0.1); acc_dt += 0.1;
306 bal.getPosition(pos);
314 if ((alt < 2990) && (hysteresis == true))
319 if ((bal.getTemperature()-273.16)>250.0)
320 burner = false; //emergency
330 //printf("Position: (%f/%f/%f), dP: (%f/%f/%f), burner: ", pos[0], pos[1], pos[2], dp[0], dp[1], dp[2]);
331 printf("%f \t%f \t%f \t%f\n", acc_dt/60.0, bal.getTemperature()-273.16, pos[2], bal.getFuelLeft());
335 bal.set_burner_strength(1.0);
340 bal.set_burner_strength(0.0);