Моделирование полёта космической капсулы в заданную точку
Автор
sergey100779Краткая аннотация. Представлен проект моделирования способов формирования траекторий спускаемого аэробаллистического летательного аппарата (ЛА) на примере спускаемой космической капсулы требуемых конфигураций при наведении в заданную точку земной поверхности, включающий задание требуемого направления вектора конечной скорости ЛА в точке приземления, построение целевой системы координат, связанной с заданным направлением конечной скорости, формирование программ требуемых ускорений ЛА в проекциях на картинную плоскость наведения, нормальной к направлению конечной скорости, нахождение закона изменения параметров управления спускаемым ЛА из динамических уравнений движения ЛА с помощью бортовой навигационно-измерительной системы. Проект позволяет производить оценку вектора пространства состояния, динамики изменения перегрузок, скоростей, аэродинамических характеристик ЛА с целью подбора соотвтествующей траектории спуска на Землю либо наоборот,- ухода и уклонения от средств перехвата. Как вариант, проект может быть преобразован и адаптирован для моделирования полета ракетного комплекса "Орешник" на конечном участке троектории подлёта к цели при представлении разработчиком аэродинамических характеристик и требований к конечной боевой ступени.
Концепция синтеза законов управления подвижными объектами как решение обратной задачи динамики является весьма общей и подробно рассмотрена в [1-10]. Применение этой концепции к подвижным объектам приводит к методу управления, получившему название «метод требуемых ускорений» [7]. Общая структура и формульные схемы алгоритмов управления аэробаллистических летательных аппаратов (ЛА) по методу требуемых ускорений описаны в [8]. Основу метода требуемых ускорений составляет идея задания желаемого (требуемого) закона движения управляемого объекта в виде так называемых “программ требуемых ускорений”, сформированных по краевым условиям управления и критериям оптимальности. Для полностью управляемых объектов необходимо задание трех программ требуемых ускорений (по числу степеней свободы поступательного движения объекта) в проекции на оси той или иной прямоугольной системы координат. Для неполностью управляемых объектов (таких, как аэробаллистические ЛА) достаточно задания двух программ требуемых ускорений. В частности, в задачах наведения аэробаллистических ЛА в заданную точку прицеливания программы требуемых ускорений задаются в проекциях на две оси целевой системы координат, лежащих в так называемой картинной плоскости, нормальной к направлению вектора терминальной скорости аэробаллистического ЛА в точке прицеливания [10]. Для завершения синтеза управлений по методу требуемых ускорений осуществляется приравнивание программных ускорений правым частям динамических уравнений движения ЛА и из полученных уравнений находятся искомые значения параметров управления объекта, обеспечивающие движения по траектории, заданной программами требуемых ускорений. Отметим важнейшие особенности метода требуемых ускорений:
- Метод позволяет осуществлять синтез законов терминального управления подвижными объектами при нелинейных моделях движения, не требуя линеаризации или каких-либо других упрощений этих моделей.
- Качественные характеристики синтезированных законов управления (оптимальность, точность управления, оперативность подготовки данных на пуск) определяются главным образом видом применяемых программ требуемых ускорений. В частности, замкнутость программ ускорений (т.е. задание их в виде функций текущих параметров движения объекта) определяет свойство замкнутости законов управления в целом, чем обеспечивается компенсация внешних и внутренних возмущающих воздействий на объект и высокая точность управления.
In [ ]:
#Ввод исходных данных
C_x=0.1#Аэродинамические коэффициенты упрощенная модель константы
C_y=0.3
C_z=0.3
S=2 #Площадь поперечного сечения спасательной капсулы
m=500 #Масса
k1=1 #Коэффициент для программы управления
k2=1 #Коэффициент для программы управления
R_zem=6371000 #Радиус Земли
ro_0=1.225
Pi_0=3.986e+14
h_m=7110
dt=0.5 #Интрвал времени расчета
V_n0=6500 #Начальная скорость в м/с
h_n=50000 #Высота входа в атмосферу в м
#Географические широта и долгота 1град = 111194 м = 111 км
fi_n=45
fi_z=45.6182964482
ljambda_n=0
ljambda_z=0
tetta_n=-23 #Угол бросания
tetta_z=-16 #Угол подлёта к цели
#Угол рыскания
psi_n=0
psi_z=0
#Преобразование значений углов
fi_n=pi*fi_n/180
fi_z=fi_z*pi/180
ljambda_n=ljambda_n*pi/180
ljambda_z=ljambda_z*pi/180
tetta_n=tetta_n*pi/180
tetta_z=tetta_z*pi/180
psi_n=psi_n*pi/180
psi_z=psi_z*pi/180
e_n=[sin(tetta_z)*cos(psi_z)*sin(fi_z)*sin(ljambda_z)+cos(tetta_z)*cos(fi_z)*sin(ljambda_z)+sin(tetta_z)*sin(psi_z)*cos(ljambda_z)
-sin(tetta_z)*cos(psi_z)*cos(fi_z)+cos(tetta_z)*sin(fi_z)
sin(tetta_z)*cos(psi_z)*sin(fi_z)*cos(ljambda_z)+cos(tetta_z)*cos(fi_z)*cos(ljambda_z)-sin(tetta_z)*sin(psi_z)*sin(ljambda_z)]
e_b=[-sin(psi_z)*sin(fi_z)*sin(ljambda_z)+cos(psi_z)*cos(ljambda_z)
sin(psi_z)*cos(fi_z)
-sin(psi_z)*sin(fi_z)*cos(ljambda_z)-cos(psi_z)*sin(ljambda_z)]
r_z=R_zem*[cos(fi_z)*sin(ljambda_z) sin(fi_z) cos(fi_z)*cos(ljambda_z)]
#Ввод начальных значений
r=R_zem+h_n
fi=fi_n
ljambda=ljambda_n
V=V_n0
tetta=tetta_n
psi=psi_n
t=0
i=1
r_0=0
fi_0=0
ljambda_0=0
V_0=0
tetta_0=0
psi_0=0
# Матрица
A_z_T=[cos(tetta_z)*cos(psi_z) sin(tetta_z) -cos(tetta_z)*sin(psi_z);
-sin(tetta_z)*cos(psi_z) cos(tetta_z) sin(tetta_z)*sin(psi_z);
sin(psi_z) 0 cos(psi_z)]
# Матрица
A_T_G=[-sin(ljambda_z)*sin(fi_z) cos(fi_z) -cos(ljambda_z)*sin(fi_z);
sin(ljambda_z)*cos(fi_z) sin(fi_z) cos(ljambda_z)*cos(fi_z);
cos(ljambda_z) 0 -sin(ljambda_z)]
#Начальный блок вывода параметров
ro=ro_0*exp(-(r-R_zem)/h_m)
q=(ro*V^2)/2
g=Pi_0/(r^2)
r_t=r*[cos(fi)*sin(ljambda) sin(fi) cos(fi)*cos(ljambda)]
delta_r=r_t-r_z
ro_n=delta_r*e_n
ro_b=delta_r*e_b
#Вектор контролируемых параметров с переходом в массив WR
WR=zeros(32000,14)
WR[i,1]=t
WR[i,2]=V
WR[i,3]=0 #n_x
WR[i,4]=0 #n_y
WR[i,5]=180*tetta/pi
WR[i,6]=0
WR[i,7]=180*fi/pi
WR[i,8]=r-R_zem
WR[i,9]=(fi_z-fi)*R_zem
WR[i,10]=q
WR[i,11]=0 #ro_n
WR[i,12]=0 #ro_b
WR[i,13]=0 #ro_v_j зарезервирован
WR[i,14]=0 #ro_v_k зарезервирован
#Конец начального блока вывода параметров
#Рабочий цикл интегрирования параметров управления полётом
flag=1
while flag==1
x=(fi-fi_n)*R_zem
y=(ljambda-ljambda_n)*R_zem
ro=ro_0*exp(-(r-R_zem)/h_m)
q=(ro*V^2)/2
g=Pi_0/(r^2)
V_t=V*[-cos(tetta)*cos(psi)*sin(fi)*sin(ljambda)+sin(tetta)*cos(fi)*sin(ljambda)-cos(tetta)*sin(psi)*cos(ljambda) cos(tetta)*cos(psi)*cos(fi)+sin(tetta)*sin(fi) -cos(tetta)*cos(psi)*sin(fi)*cos(ljambda)+sin(tetta)*cos(fi)*cos(ljambda)+cos(tetta)*sin(psi)*sin(ljambda)]
r_t=r*[cos(fi)*sin(ljambda) sin(fi) cos(fi)*cos(ljambda)]
delta_r=r_t-r_z
ro_n=delta_r*e_n
ro_b=delta_r*e_b
V_n=V_t*e_n
V_b=V_t*e_b
V_t_delta_r=V_t*transpose(delta_r)
T_=(norm(delta_r)^4)/abs(V_t_delta_r[1])
a_tr_n=-((2+3*k1+k1^2)*ro_n)/((T_)^2)-(2*(k1+1)*V_n)/(T_)#Программа требуемых ускорений
a_tr_b=-((2+3*k2+k2^2)*ro_b)/((T_)^2)-(2*(k2+1)*V_b)/(T_)#Программа требуемых ускорений
# Матрица
A_G_t=[-sin(ljambda)*sin(fi) sin(ljambda)*cos(fi) cos(ljambda);
cos(fi) sin(fi) 0;
-cos(ljambda)*sin(fi) cos(ljambda)*cos(fi) -sin(ljambda)]
# Матрица
A_t_ps=[cos(tetta)*cos(psi) -sin(tetta)*cos(psi) sin(psi);
sin(tetta) cos(tetta) 0;
-cos(tetta)*sin(psi) sin(tetta)*sin(psi) cos(psi)]
A_z_ps=A_z_T*A_T_G*A_G_t*A_t_ps
Q_=C_x*S*q
f_1=a_tr_n[1]+(Q_/m)*A_z_ps[2,1]+g*cos(tetta_z)
f_2=a_tr_b[1]+(Q_/m)*A_z_ps[3,1]
D=A_z_ps[2,2]*A_z_ps[3,3]-A_z_ps[2,3]*A_z_ps[3,2]
alpha=(m/(C_y*q*S))*((f_1*A_z_ps[3,3]-f_2*A_z_ps[2,3])/D)
if alpha>=pi/6
alpha=pi/6
else alpha=alpha
end
if alpha<=-pi/6
alpha=-pi/6
else alpha=alpha
end
betta=(m/(C_z*q*S))*((f_1*A_z_ps[3,2]-f_2*A_z_ps[2,2])/D)
if betta>=pi/6
betta=pi/6
else betta=betta
end
if betta<=-pi/6
betta=-pi/6
else betta=betta
end
X=-C_x*S*q
Y=C_y*alpha*S*q
Z=-C_z*betta*S*q
n_x=(C_x*S/m)*(q/g)
n_y=Y/(g*m)
n_z=Z/(g*m)
# Решение системы дифференциальных уравнений Рунге-Кутта 4 порядка
r_0=r
fi_0=fi
ljambda_0=ljambda
V_0=V
tetta_0=tetta
psi_0=psi
X_0=[r_0 fi_0 ljambda_0 V_0 tetta_0 psi_0]
f_r=V*sin(tetta)
f_fi=(V/r)*cos(psi)*cos(tetta)
f_ljambda=-(V/r)*sin(psi)*cos(tetta)/cos(fi)
f_V=-(C_x*S*ro_0*exp(-(r-R_zem)/h_m)*(V^2/2))/m-(Pi_0/(r^2))*sin(tetta)
f_tetta=(alpha*C_y*S*ro_0*exp(-(r-R_zem)/h_m)*(V^2/2))/(m*V)-Pi_0/(V*(r^2))*cos(tetta)+(V/r)*cos(tetta)
f_psi=(betta*C_z*S*ro_0*exp(-(r-R_zem)/h_m)*((V^2)/2))/(m*V*cos(tetta))+(V/r)*tan(fi)*sin(psi)*cos(tetta)
k_1=[f_r*dt f_fi*dt f_ljambda*dt f_V*dt f_tetta*dt f_psi*dt]
X_=X_0+k_1/2
#--------------
r=X_[1,1]
fi=X_[1,2]
ljambda=X_[1,3]
V=X_[1,4]
tetta=X_[1,5]
psi=X_[1,6]
f_r=V*sin(tetta)
f_fi=(V/r)*cos(psi)*cos(tetta)
f_ljambda=-(V/r)*sin(psi)*cos(tetta)/cos(fi)
f_V=-(C_x*S*ro_0*exp(-(r-R_zem)/h_m)*(V^2/2))/m-(Pi_0/(r^2))*sin(tetta)
f_tetta=(alpha*C_y*S*ro_0*exp(-(r-R_zem)/h_m)*(V^2/2))/(m*V)-Pi_0/(V*(r^2))*cos(tetta)+(V/r)*cos(tetta)
f_psi=(betta*C_z*S*ro_0*exp(-(r-R_zem)/h_m)*((V^2)/2))/(m*V*cos(tetta))+(V/r)*tan(fi)*sin(psi)*cos(tetta)
k_2=[f_r*dt f_fi*dt f_ljambda*dt f_V*dt f_tetta*dt f_psi*dt]
X_=X_0+k_2/2
#--------------
r=X_[1,1]
fi=X_[1,2]
ljambda=X_[1,3]
V=X_[1,4]
tetta=X_[1,5]
psi=X_[1,6]
f_r=V*sin(tetta)
f_fi=(V/r)*cos(psi)*cos(tetta)
f_ljambda=-(V/r)*sin(psi)*cos(tetta)/cos(fi)
f_V=-(C_x*S*ro_0*exp(-(r-R_zem)/h_m)*(V^2/2))/m-(Pi_0/(r^2))*sin(tetta)
f_tetta=(alpha*C_y*S*ro_0*exp(-(r-R_zem)/h_m)*(V^2/2))/(m*V)-Pi_0/(V*(r^2))*cos(tetta)+(V/r)*cos(tetta)
f_psi=(betta*C_z*S*ro_0*exp(-(r-R_zem)/h_m)*((V^2)/2))/(m*V*cos(tetta))+(V/r)*tan(fi)*sin(psi)*cos(tetta)
k_3=[f_r*dt f_fi*dt f_ljambda*dt f_V*dt f_tetta*dt f_psi*dt]
X_=X_0+k_3/2
#---------------
r=X_[1,1]
fi=X_[1,2]
ljambda=X_[1,3]
V=X_[1,4]
tetta=X_[1,5]
psi=X_[1,6]
f_r=V*sin(tetta)
f_fi=(V/r)*cos(psi)*cos(tetta)
f_ljambda=-(V/r)*sin(psi)*cos(tetta)/cos(fi)
f_V=-(C_x*S*ro_0*exp(-(r-R_zem)/h_m)*(V^2/2))/m-(Pi_0/(r^2))*sin(tetta)
f_tetta=(alpha*C_y*S*ro_0*exp(-(r-R_zem)/h_m)*(V^2/2))/(m*V)-Pi_0/(V*(r^2))*cos(tetta)+(V/r)*cos(tetta)
f_psi=(betta*C_z*S*ro_0*exp(-(r-R_zem)/h_m)*((V^2)/2))/(m*V*cos(tetta))+(V/r)*tan(fi)*sin(psi)*cos(tetta)
k_4=[f_r*dt f_fi*dt f_ljambda*dt f_V*dt f_tetta*dt f_psi*dt]
d_X_=(k_1+2*k_2+2*k_3+k_4)/6
X_=X_0+d_X_
r=X_[1,1]
fi=X_[1,2]
ljambda=X_[1,3]
V=X_[1,4]
tetta=X_[1,5]
psi=X_[1,6]
r=r_0+d_X_[1,1]
fi=fi_0+d_X_[1,2]
ljambda=ljambda_0+d_X_[1,3]
V=V_0+d_X_[1,4]
tetta=tetta_0+d_X_[1,5]
psi=psi_0+d_X_[1,6]
r_0=r
fi_0=fi
ljambda_0=ljambda
V_0=V
tetta_0=tetta
psi_0=psi
i=i+1
WR[i,1]=t
WR[i,2]=V
WR[i,3]=n_x
WR[i,4]=n_y
WR[i,5]=180*tetta/pi
WR[i,6]=180*alpha/pi
WR[i,7]=180*fi/pi
WR[i,8]=r-R_zem
WR[i,9]=(fi_z-fi)*R_zem
WR[i,10]=q
WR[i,11]=ro_n[1]
WR[i,12]=V_n[1]
WR[i,13]=T_
WR[i,14]=0 #a_tr_n[1]/9.8
t=t+dt
if R_zem>=r
flag=0
n=i
end
end
#конец основного цикла
#Построение графиков
x = WR[1:i-1,1]
y1 = @.WR[1:i-1,2]
y2 = @.WR[1:i-1,8]
y3 = @.WR[1:i-1,4]
y4 = @.WR[1:i-1,10]
y = [y1 y2 y3 y4]
plot(x,y, layout=(4, 1), label=["V-скорость от времени" "H-высота от времени" "n_y-перегрузка от времени" "q-скоростной напор от времени"], legend=true)
Out[0]:
Выводы
Представленный проект позволяет производить анализ возможных подходов к видоизменению и адаптации алгоритмов управления движением ЛА на конечном участке траектории на основе метода требуемых ускорений. Данные алгоритмы объединяют в себе одновременно две разнородные задачи маневрирования и наведения.
- Алгоритмы по методу требуемых ускорений содержат в себе обширный ряд свободных параметров, позволяющих случайным образом заранее предусматривать и задавать широкий спектр напряженных траекторий маневров. Отличительной особенностью алгоритмов наведения по методу требуемых ускорений является необходимость априорного задания требуемого направления вектора конечной скорости ЛА в точке цели, при этом к величине скорости вследствие неуправляемости ЛА по параметрам продольного движения требования не предъявляются. По заданному направлению конечного вектора скорости определяется ориентация нормальной к ней картинной плоскости наведения, в которой задаются программы требуемых ускорений. Движение аэробаллистического ЛА в проекции на картинную плоскость полностью управляемо с помощью подъемной и боковой аэродинамических сил. Направление вектора конечной скорости ЛА задается двумя углами – углом бросания и курсовым углом. Данные углы могут рассматриваться в качестве свободных параметров, выбором которых в допустимых пределах, определяемых характеристиками маневренности ЛА, возможно формирование различных траекторий маневра. Так, выбором угла бросания при неизменных начальных условиях движения ЛА возможно формирование траекторий маневра типа “кабрирование-пикирование”, “пикирование-кабрирование” и изменение знака поперечной перегрузки аппарата в процессе маневра. Выбором курсового угла возможно формирование траекторий бокового маневра ЛА (по отношению к плоскости баллистического спуска), обеспечивающих подлет к цели с заданного направления. В общем случае углы могут изменяться в процессе наведения по заданной программе. В частности, программное изменение углов, при котором вектор конечной скорости вращается по образующей кругового конуса с вершиной в точке цели, обеспечивает формирование спиралевидных траекторий маневра ЛА при наведении [6]. В данном случае в качестве свободных параметров управления маневром могут рассматриваться угол полураствора упомянутого конуса и период вращения вектора конечной скорости по образующей этого конуса. Другую группу свободных параметров в алгоритмах наведения образуют константы в программах требуемых ускорений. В качестве третьей группы свободных параметров в алгоритмах наведения ЛА могут рассматриваться весовые множители в критерии оптимальности, которыми определяются вид программ требуемых ускорений. Изменение вида применяемых программ ускорений очевидным образом влияет на форму траекторий маневра ЛА при наведении.
- Синтезированы пять вариантов программ управления требуемыми ускорениями (в проекте не представлены) алгоритмов по методу требуемых ускорений в задаче с обобщенным квадратичным критерием по принципу максимума Л.С. Понтрягина. Данные варианты программ требуемых ускорений содержат в себе экспоненциальные и периодические функции, позволяя, тем самым, формировать управление для выполнения сложнопрогнозируемых траектории различного профиля.
- При помощи проекта получены способы формирования траекторий противоракетного маневра ЛА требуемого профиля посредствам изменения вектора терминальной скорости. Изменение вектор скорости может задаваться как постоянным с целью формирования траекторий кабрирования, пикирования, подлетных с заданного направления, так и программным образом в виде вращения вектора скорости по образующей конуса.
- На основе метода требуемых ускорений разработан программно-терминальный вариант алгоритма наведения ЛА на конечном участке траектории. Данный вариант позволяет существенно расширить совокупность возможных специфических траекторий с учетом наложенных ограничений на объект управления (в частности для ЛА планирующего типа). Первоначальное моделирование движения ЛА показало, что методическая погрешность наведения при шаге интегрирования по времени до 1 се¬кунд и выполнению профиля программно-терминальной траектории являлась не выше 1 метра. Надлежащим выбором шага интегрирования и подбором опорных точек достигается ме¬тодическая погрешность наведения менее 1 м.
- Весьма полезны способы модернизации и адаптации алгоритмов по методу требуемых ускорений для выполнения специфических задач. Введение коэффициента адаптации и коррекция времени прогноза в программы требуемых ускорений позволяет избавиться от высокого уровня значений программ требуемых ускорений, выводящих органы управления объекта на придельные значения и лишающих возможности совершать маневры сложного профиля.
Список литературы
- Александер Т.М. Наведение управляемого летательного аппарата на фиксированную на Земле целевую точку.- NSWS, TR-3916, USA, 1978.
- Батенко А.П. Системы терминального управления. М.: Радио и связь, 1984. 160 с.
- Батенко А.П. Управление конечным состоянием движущихся объектов. М.: Сов. Радио, 1977.256 с.
- Бородовский В. Н. Управление конечными параметрами движения летательных аппаратов. М.: МО РФ, 2005. 232 с.
- Оссовский С., Нейронные сети для обработки информации. М.: «Финансы и статистика», 2004.
- Патент 2296940 C1 РФ, МПК F41G 7/22, F41G 7/34 Способ формирования траекторий спускаемого аэробаллистического ЛА требуемых конфигураций при наведении в заданную точку земной поверхности / Г.Н. Разорёнов, С.Н. Шевцов (РФ). 2005134996/28; Заяв. 11.11.2005; Опубл. 10.04.07 Бюл. №10.
- Петров Б.Н., Портнов-Соколов Ю.П., Андриенко А.Я., Иванов В.П. Бортовые терминальные системы управления. М.: Машиностроение, 1983.
- Разоренов Г.Н., Бахрамов Э.А., Титов Ю.Ф. Системы управления летательными аппаратами. М.: Машиностроение, 2003.
- Разорёнов Г.Н., Самарин А.А. Теория и системы оптимального управления М.: ВА РВСН, 2007.
- Разорёнов Г.Н., Шевцов С.Н. Алгоритмы управления движением аэробаллистических ЛА, обеспечивающие высокоточное наведение на цель и формирование траекторий заданного профиля / Труды МИТ Т. 8/ М.: МИТ, 2006. с. 102-128.