Реферат: Математическое моделирование космических систем

4. Задача. Математическое моделирование движения ракеты-носителя

Задание:

Разработать программу и провести расчеты в среде MatLab + Simulink по исследованию математической модели движения ракеты-носителя в плотных слоях атмосферы при выведении полезного груза на орбиту. В качестве исходных данных для моделирования принять материалы из Википедии (Интернет) по ракете-носителю "Союз-2".

Решение задания:

В среде MatLab + Simulink разработана программа, схема которой (mdl- file) представлена на рисунке 1.

Рисунок 1 - Схема программы моделирования движения РН в плотных слоях атмосфере (mdl- file)

Математическое моделирование движения ЛА основывается на интегрировании системы дифференциальных уравнений, описывающих законы изменения основных параметров ЛА при движении в плотных слоях атмосферы в системах координат OXYZ, O1X1Y1Z1 и O1XvYvZv. OXYZ - начальная стартовая система координат с началом координат O → на поверхности Земли в точке (j0,l0), где j0 и l0 - широта и долгота точки старта, плоскость OXZ → в плоскости местного горизонта, ось OX → по направлению расчетной точки выведения (точки падения, взлета и т.п.) ЛА, ось OY → по направлению местной вертикали. O1X1Y1Z1 - связанная система координат с началом координат O1 → в центре масс ЛА, ось O1X1 → по направлению продольной оси ЛА, ось O1Y1 → в поперечной плоскости по направлению оси симметрии аппарата. O1XvYvZv - скоростная система координат с осью OvXv → по направлению воздушного потока. Угловая ориентация систем координат (OXYZ) и (O1X1Y1Z1) определяется углами тангажа - u, рысканья - y и крена - g. Угловая ориентация систем координат (O1X1Y1Z1) и (O1XvYvZv) определяется углами атаки - a, скольжения - b и аэродинамического угла крена - f.

Переход от начальной стартовой системы координат (OXYZ) к связанной системе координат (O1X1Y1Z1) может быть осуществлён с помощью матрицы:

следующим образом: y1= А·y, где ψ, υ, γ - углы рысканья, тангажа и крена. А=с(γ)·в(υ)·а(ψ), где с,в,а матрицы последовательных поворотов вектора y на углы ψ, υ, γ. y и y1 - вектора положения в системах координат OXYZ и O1X1Y1Z1, соответственно.

|cos(ψ)0-sin(ψ)|

а = |010 |

|sin(ψ)0 cos(ψ)| ,

| cos(υ)sin(υ)0|= |- sin(υ)cos(υ)0|

| 100 |

c = | 1cos(γ)sin(γ) |

|0 - sin(γ) cos(γ) | .

Уравнения движения центра масс ЛА представлены в следующем виде:

m× ÿ = AТ×(T + Av×Ra) + G,

орбита спутник математический моделирование

где ÿ - вектор ускорения ЛА y={x,y,z}, T - вектор тяги маршевых двигателей ЛА, Ra={Xa,Ya,Za} - вектор аэродинамических сил в поточной с.к.( Xa=Сxa·Sm·q, Ya=Сya·Sm·q, Za=Сza·Sm·q, q=ρV2/2 - скоростной напор, ρ - плотность, Sm - площадь миделевого сечения, Сxa, Сya, Сza - аэродинамические коэффициенты, V - модуль вектора воздушной скорости), G - вектор силы тяготения, AТ - матрица перехода от связанной системы координат (с.к.) к начальной стартовой с.к., определяемая углами u,y и g, Av - матрица перехода от поточной с.к. к связанной с.к., определяемая углами a и b, m = m0 - dm× t (m0 - начальная масса ЛА, dm - массовый секундный расход, t - текущее время полета). Вектор воздушной скорости ЛА, используемый при определении углов a и b, вычисляется следующим образом Vv =V + W. V - вектор скорости ЛА, определяемый, как и текущий вектор положения ЛА - r, в результате интегрирования дифференциальных уравнений линейного движения ЛА, W - текущий вектор скорости ветра.

Вектор угловой скорости объекта можно представить в виде суммы следующих слагаемых: ω. Проектируя это уравнение на оси связанной системы координат, получим кинематические соотношения:

Уравнения вращательного движения:

J·ώ=Ma +Mu,

где Mu и Ma моменты сил управления и аэродинамических сил, определяемых в т.ч. текущими отклонениями центров масс и давления ЛА. В расчётах принята экспоненциальная модель атмосферы: r = r0×е-b×h, где b - градиент изменения плотности по высоте полёта ЛА- h.

Для заданной программы движения ЛА в плотных слоях атмосферы на рисунке 2 представлен график изменения тяги двигателей на участке выведения.

Рисунок 2 - График изменения тяги двигателей на участке выведения

Задание 2:

Доработать, программу моделирования движения ракеты-носителя программным блоком, позволяющим определять параметры орбиты выведения полезной нагрузки. Представить в табличном или графическом виде данные параметрам орбиты выведения в зависимости от начальных условий старта (по долготе и широте точки старта, азимуту плоскости выведения) и времени выключения двигателей последней ступени ракеты-носителя.

Решение задания:

Параметры орбиты выведения определяются терминальными параметрами движения ракеты-носителя, определяемыми следующими соотношениями:

(Представить основные уравнения эллиптического движения и зависимости элементов орбиты с начальными условиями - терминальными (конечными) параметрами выведения ракеты-носителя).

В среде MatLab + Simulink разработана блок программы, m- file которой представлен на рисунке 3.

function [f]=fnomin2(q);=6371000;muZ=3.98602*10^14;=0.0;fi=1.0;az=0.0;(1,1)=q(1);q0(2,1)=q(2)+Rz;q0(3,1)=q(3);v(1,1)=q(4);q0v(2,1)=q(5);q0v(3,1)=q(6);(1,1)=0;u(1,2)=0;u(1,3)=-1;(2,1)=0;u(2,2)=1;u(2,3)=0;(3,1)=1;u(3,2)=0;u(3,3)=0;(1,1)=1;a(1,2)=0;a(1,3)=0;(2,1)=0;a(2,2)=sin(1);a(2,3)=cos(1);(3,1)=0;a(3,2)=-cos(1);a(3,3)=sin(1);(1,1)=cos(fi);b(1,2)=sin(fi);b(1,3)=0;(2,1)=-sin(fi);b(2,2)=cos(fi);b(2,3)=0;(3,1)=0;b(3,2)=0;b(3,3)=1;(1,1)=cos(az);c(1,2)=0;c(1,3)=-sin(az);(2,1)=0;c(2,2)=1;c(2,3)=0;(3,1)=sin(az);c(3,2)=0;c(3,3)=cos(az);=u*a*b*c;qg=uabc*q0; qgv=uabc*q0v;=sqrt(qg(1,1)^2+qg(2,1)^2+qg(3,1)^2);vg=sqrt(qgv(1,1)^2+qgv(2,1)^2+qgv(3,1)^2);vkr=sqrt(muZ/rg);=acos((qg(1,1)*qgv(1,1)+qg(2,1)*qgv(2,1)+qg(3,1)*qgv(3,1))/rg/vg);=qg(2,1)*qgv(3,1)-qg(3,1)*qgv(2,1);c2=qg(3,1)*qgv(1,1)-qg(1,1)*qgv(3,1);c3=qg(1,1)*qgv(2,1)-qg(2,1)*qgv(1,1);=-muZ*qg(1,1)/rg+c3*qgv(2,1)-c2*qgv(3,1);f2=-muZ*qg(2,1)/rg+c1*qgv(3,1)-c3*qgv(1,1);f3=-muZ*qg(3,1)/rg+c2*qgv(1,1)-c1*qgv(2,1);=sqrt(c1^2+c2^2+c3^2);fs=sqrt(f1^2+f2^2+f3^2);=cs^2/muZ;e=fs/muZ;ra=p;rp=p;(e<1) ra=p/(1-e);rp=p/(1+e);end;=(rp+ra)/2;=2*pi*sqrt(ap^3/muZ);=acos(c3/cs);=0;if in>0 kom=-c2/(cs*sin(in));som=c1/(cs*sin(in));om=acos(kom);end;=qg(1,1)/rg;mr=qg(2,1)/rg;nr=qg(3,1)/rg;lv=qgv(1,1)/vg;mv=qgv(2,1)/vg;nv=qgv(3,1)/vg;(1)=ra;(2)=rp;(3)=p;(4)=e;(5)=in;(6)=om;(7)=vkr;

Таблица 1 - результатов моделирования параметров орбиты выведения

Список использованных источников

1.  Ю.Г. Сихарулидзе Баллистика летательных аппаратов. Москва, "Наука", 1982.

2.      Л.Н. Лысенко Наведение и навигация баллистических ракет. Издательство МГТУ им. Н.Э. Баумана, 2007 г.

.        А.П. Разыграев Основы управления полётом космических аппаратов и кораблей.