Введение
Определение скорости подъёма и спуска летательных аппаратов легче воздуха (ЛАЛВ) до настоящего времени является практически важной задачей, возникающей при проектировании таких аппаратов.
Большое количество публикаций посвящено ЛАЛВ, например, только на нашем ресурсе приведены две очень интересные статьи [1,2], касающиеся истории развития на примере конкретных конструкций дирижаблей и стратостатов. Однако очень мало расчётов динамики вертикального полёта таких устройств, позволяющих хотя бы ориентировочно определять скорости подъёма и спуска ЛАЛВ.
Последнее утверждение требует определённого пояснения, поскольку искушённый читатель хорошо помнит школьный курс физики, в котором решались задачи на высоту подъёма и другие параметры воздушных шаров, заполненных газами легче воздуха или самим подогреваемым во время полёта воздухом.
Все указанные задачи были основаны на равенстве двух сил: силы веса и выталкивающей силы. Газы считались идеальными и их параметры вычислялись по закону Менделеева Клапейрона. Однако, даже простой учёт третьей силы сопротивления воздуха уже приводит к системе дифференциальных уравнений, которая аналитически не решается. Необходимо так же учитывать изменение плотности атмосферного воздуха с высотой подъёма и температурой.
Кроме этого, если нужно рассмотреть не только подъём, но и зависание шара и его спуск на землю, то совсем уж не детская задача получается. Надеюсь, что рассмотрение решения подобной задачи средствами Python не только будет способствовать расширению знаний по физике, но и популяризации самого языка программирования Python. Что я и пытаюсь делать в своих публикациях на этом ресурсе.
Математическая модель полёта ЛАЛВ с оболочкой в форме шара, объём которого не изменяется с изменением высоты
Ограничимся рассмотрением движения его центра масс под действием следующих сил: силы тяжести (G), архимедовой силы (Fa) и силой аэродинамического сопротивления (Fc). Запишем соотношения для определения сил через параметры движения и воздушной среды[3]:
В приведенных формулах приняты обозначения: h – высота подъема шара, dh/dt – вертикальная скорость, m – масса, g – ускорение свободного падения, W – объем шара, c – коэффициент лобового сопротивления, S – характерная площадь сопротивления (площадь миделя).
Зависимость плотности воздуха от высоты будем полагать экспоненциальной:
где
– плотность воздуха на нулевой высоте, b–коэффициент. Сила тяжести направлена вниз, архимедова сила – вверх, а сила аэродинамического сопротивления всегда направлена «против движения», поэтому корректный учет этой силы в уравнениях движения требует введения множителя 
.Однако, для наших целей этот факт не имеет принципиального значения, и мы ограничимся рассмотрением только этапа подъема шара, когда сила аэродинамического сопротивления направлена вниз и, следовательно, будет учтена в уравнениях движения со знаком минус. Теперь уравнение движения может быть записано в виде:
, (1)Дополнительно предположим, что воздушный шар представляет собой однородное тело радиуса R с плотностью
. Тогда величина площади, определяющая его аэродинамическое сопротивление, определится как 
, объем как 
, а масса, соответственно, как 
.Теперь видно, что каждый член уравнения (1) содержит в качестве множителя величину S. Следовательно, каждый член уравнения движения может быть сокращен на величину множителя S, а само уравнение примет вид:
, (2)Введём обозначения:
;
;и перепишем (2) в виде следующей системы нелинейных уравнений:
, (3)Влияние на скорость и высоту подъёма ЛАЛВ температуры атмосферного воздуха
Для этого сначала решим систему (3) с использованием следующего соотношения для зависимости плотности атмосферного воздуха от высоты без учёта температуры:
Повторим решение системы (3), но уже с использованием соотношения для зависимости плотности воздуха от высоты и температуры:
где: b=0.000125 — константа, связанная с плотностью воздуха в 1/м.;
a=0.0065 — константа, связанная с температурой воздуха в K/м.
– температура на уровне моря.Листинг программы
# -*- coding: utf8 -*-
from numpy import*
from scipy.integrate import odeint
import matplotlib.pyplot as plt
g=9.81# ускорение свободного падения на земле в м/с2.
rv=1.29# плотность атмосферного воздуха в кг/м3.
rg=0.17# плотность гелия в кг/м3.
R=8# радиус оболочки ЛАЛВ в м.
b=0.000125# константа, связанная с плотностью воздуха в 1/м
a=6.5*10**-3# константа, связанная с температурой воздуха в К/м
c=0.4#коэффициент лобового сопротивления
mo=240#масса в кг
V=(4/3)*pi*R**3
rs=rg+mo/V# суммарная плотность материала ЛАЛВ, массы гелия, и нагрузки
p1=rv/rs# введенный параметр
p2=3*c/(8*R)# введенный параметр
T0=300
def fun(y, t):
y1, y2= y
return [y2,-g+g*p1*exp(-b*y1*T0/(T0-a*y1))-p1*p2*exp(-b*y1*T0/(T0-a*y1))*y2**2]
t =arange(0,1100,0.01)
y0 = [0.0,0.0]
[y1,y2]=odeint(fun, y0,t, full_output=False).T
plt.title("Характеристики подъёма ЛАЛВ \n Объём: %s м3. Масса : %s кг. \n Подъёмная сила: %s kН. "%(round(V,0),mo,round(0.001*g*rv*V,0)))
plt.plot(t/60,y1,label='Максимальная высота подъёма: %s км. \n Максимальная скорость: % s м/с .\n С учётом температуры воздуха'%(round(max(y1)/1000,2),round(max(y2),2)))
def fun(y, t):
y1, y2= y
return [y2,-g+g*p1*exp(-b*y1)-p1*p2*exp(-b*y1)*y2**2]
[y1,y2]=odeint(fun, y0,t, full_output=False).T
plt.plot(t/60,y1,label='Максимальная высота подъёма: %s км. \n Максимальная скорость: % s м/с \n Без учёта температуры воздуха'%(round(max(y1)/1000,2),round(max(y2),2)))
plt.ylabel('Высота в м')
plt.xlabel(' Время в минутах')
plt.legend(loc='best')
plt.grid(True)
plt.show()Получим:
Расчётное значение высоты подъёма ЛАЛВ с учётом температуры меньше, чем без учёта. Скорость подъёма аппарата при этом остаётся неизменной.
Определение характеристик всех фаз полёта ЛАЛВ от старта до приземления
Для построения программы полёта ЛАЛВ рассмотрим условия для следующих периодов времени:
Подъём —
;Зависание —
;Приземление —
.Листинг программы
# -*- coding: utf8 -*-
from numpy import*
from scipy.integrate import odeint
import matplotlib.pyplot as plt
g=9.81# ускорение свободного падения на земле в м/с2.
rv=1.29# плотность атмосферного воздуха в кг/м3.
rg=0.17# плотность гелия в кг/м3.
R=8# радиус оболочки стратостата в м.
b=0.000125# константа, связанная с плотностью воздуха в 1/м
a=6.5*10**-3# константа, связанная с температурой воздуха в К/м
c=0.4# коэффициент лобового сопротивления
mo=240# масса в кг
V=(4/3)*pi*R**3
p2=3*c/(8*R)# введенный параметр
T0=300# температура на уровне моря
tz=4000# время зависания в секундах
rgu=1.2# плотность образовавшейся газовой смеси после стравливания гелия в кг/м3
tz=4000# время зависания
def fun(y, t):
y1,y2= y
if y2<=0:
if t<tz:
return [y2,-g+g*(rv/(rg+mo/V))*exp(-b*y1*T0/(T0-a*y1))+(rv/(rg+mo/V))*p2*exp(-b*y1*T0/(T0-a*y1))*y2**2]
elif t>=tz:
return [y2,-g+g*(rv/(rgu+mo/V))*exp(-b*y1*T0/(T0-a*y1))+(rv/(rgu+mo/V))*p2*exp(-b*y1*T0/(T0-a*y1))*y2**2]
else:
return [y2,-g+g*(rv/(rg+mo/V))*exp(-b*y1*T0/(T0-a*y1))-(rv/(rg+mo/V))*p2*exp(-b*y1*T0/(T0-a*y1))*y2**2]
t =arange(0,tz+555,0.1)
y0 = [0.0,0.0]
[y1,y2]=odeint(fun, y0,t, full_output=False).T
plt.title("Подъём, зависание, спуск ЛАЛВ \n с жёсткой оболочкой сферической формы \n Объём: %s м3. Масса : %s кг. Подъёмная сила: %s kН. "%(round(V,0),mo,round(0.001*g*rv*V,0)))
plt.plot(t,y1,label='Максимальная высота подъёма: %s км. \n Максимальная скорость: % s м/с .\n Время зависания %s с.'%(round(max(y1)/1000,2), round(max(y2),2),tz-2*555))
plt.ylabel('Высота в м')
plt.xlabel(' Время в сек.')
plt.legend(loc='best')
plt.grid(True)
plt.show()Получим:
Как следует из приведенного графика и листинга программы, для проведения вычислительного эксперимента достаточно ввести необходимые исходные данные.
Математическая модель полёта ЛАЛВ с оболочкой, объём которой изменяется с изменением высоты
К подобным ЛАЛВ относятся стратостаты. Стратостат нельзя полностью надуть гелием, придав ему максимальную подъёмную силу, которая превратит форму его оболочки в шар. Такой шар на большой высоте может лопнуть под действием возросшей разности внутреннего и наружного давлений.
По указанным причинам для расчётов максимально достижимой высоты подъёма вводят два значения его объёма: минимальный Vmin и максимальный Vmax соответственно. С учётом введенных переменных и зависимости плотности воздуха от высоты соотношения для выталкивающей силы Fa и силы тяжести Fт примут вид:
, (4)
, (5)где: M — масса оболочки и оборудования стратостата;
— плотность гелия.Приравнивая соотношения (4) и (5), предполагая, что объем оболочки V является функцией от высоты подъёма ЛАЛВ, получим соотношение:
. (6)Численные значения параметров входящих в соотношение (6) приводятся в листинге для построения графика, который приводится только с указанной целью.
Листинг графика с данными
# -*- coding: utf8 -*-
from numpy import*
from scipy.integrate import odeint
import matplotlib.pyplot as plt
g=9.81# ускорение свободного падения на земле в м/с2.
rv=1.29# плотность атмосферного воздуха в кг/м3.
rg=0.17# плотность гелия в кг/м3.
Vmin=400# начальный объём шара в/м3.
b=0.000128# константа, связанная с плотностью воздуха в 1/м.
c=0.8#коэффициент лобового сопротивления
mo=40#масса в кг
rs=rg+mo/Vmin# суммарная плотность материала стратостата, массы гелия и нагрузки
p1=rv/rs# введенный параметр
h=[(10**-3)*log((rv*w)/(mo+rg*Vmin))*b**-1 for w in arange(1*10**3,1.8*10**5,1000)]
v=[(10**-3)*w for w in arange(1*10**3,1.8*10**5,1000)]
plt.title("Теоретическая зависимость высоты подъёма стратостата \n от его максимального объёма")
plt.plot(v,h)
plt.xlabel('Максимальный объём стратостата в тыс. м3')
plt.ylabel(' Максимальная высота подъёма в км.')
plt.grid(True)
plt.show()Получим:
Изменяя параметры ЛАЛВ, приведенные в листинге программы, можно получить приведенный график и выбрать требуемый максимальный объём оболочки при проектировании. Уточнение результатов проводят с использованием огромного опыта по созданию подобных аппаратов.
Выводы:
- Получены математические модели двух типов летательных аппаратов легче воздуха, которые позволяют проводить вычислительные экспериенты для оценочного определения параметров таких аппаратов в идеализированных условиях воздушной среды.
- Предложенная многоступенчатая схема численного решения системы дифференциальных уравнений позволяет получить вертикальную траекторию летательных аппаратов легче воздуха на этапах подъёма зависания и спуска.
Ссылки
- Пара слов про дирижабли
- На пути в космос. Стратостаты
- Рыжиков Ю.И. Современный Фортран. — СПб.: Корона принт, 2004. – 288 с.
