Проект компьютеризации расчета траектории неуправляемого летательного аппарата

Таблица. Исходные данные и полученные результаты:

Таблица. Дополнительные параметры

Таблица. Движение по направляющим

Таблица. Активный участок траектории ЛА

Таблица. Пассивный участок траектории ЛА.

1.Разработка компьютерной модели расчета траектории неуправляемого летательного аппарата

1.1 Постановка задач. Анализ и представление исходных данных

компьютерный летательный аппарат

В этой главе рассмотрены следующие задачи:

Разработка компьютерной модели расчета траектории неуправляемого ЛА, используя методы численного интегрирования Эйлера и Рунге-Кутта.

Расчет траектории пассивного участка полета неуправляемого ЛА с использованием параболической теории.

Траектория полета неуправляемого ЛА (снаряда РСЗО) состоит из трех характерных участков:

движение по направляющим;

активный участок;

пассивный участок.

Расчет дополнительных параметров:

1. Секундный массовый расход (q) рассчитан по формуле:

2. Тяга двигателя ® рассчитана по формуле:

R=Ue*q

3. Коэффициент лобового сопротивления рассчитан по формуле:

Cx=Cx (M)*i, где Cx=Cx (M)-коэффициент лобового сопротивления эталонного ЛА. Принимаем Cx =0.3.

4. Площадь миделя ЛА рассчитана по формуле:

5. Закон изменения массовой плотности воздуха от высоты полета принимаем в виде, где ro= 1.23- плотность воздуха на уровне моря; высота y в формуле берется в метрах.

6. Закон изменения массы принимаем в виде m= m0-q*t.

1.2 Модель динамики неуправляемого ЛА

1.2.1 Движение ЛА по направляющим

При расчете движения по направляющем принимается во внимание действие на ЛА следующих сил: Rтяга, Gвес ЛА, Fреакция направляющих, которая состоит из силы трения Fтр и нормальной составляющей N. Аэродинамическими силами при расчете движения по направляющим обычно пренебрегают.

Кроме того, в идеальном случае направляющие принимаются жесткими и недеформированными. Расчетная схема движения представлена на рисунке 1.1.

Рис 1.1.

Уравнение движения ЛА на направляющих в проекции на направление скорости запишется в виде где fкоэффициент трения скольжения, m=m0-Q*t — принимаемый в работе закон изменения массы ЛА, G=m*g. Записывая уравнение движения в нормальной форме Коши, получаем Последовательно интегрируя уравнение с учетом соотношения ds/dt=v, находим Полагая в формулах t=td, получаем выражение для дульной скорости Vd и длины направляющих:

Если задана длина направляющих Sh, то, решая численно уравнение, определим td. Для этого уравнение запишем в виде

Если пренебречь изменением массы при движении ЛА по направляющим, то, заменяя массу ЛА средним значением m= mср,

получаем По этим формулам можно приближенно найти vd и td.

Примечание: 1. Коэффициент трения скольжения ЛА об направляющие f равен 0.15 (f=0.15).

2. Изменением массы при движении ЛА по направляющим пренебрегаем.

1.2.2 Активный участок полета ЛА

При расчете активного участка полета ЛА принимают следующие допущения:

1. угол атаки равен нулю.

2. тяга R постоянна.

При сделанных допущениях уравнения для расчета траектории ЛА на активном участке полета запишутся в виде следующей системы дифференциальных уравнений первого порядка:

Здесь X=0.5*Cx*Ro*S*V — сила лобового сопротивления. Остальные элементы системы были указаны выше.

1.2.3 Пассивный участок полета ЛА

При расчете пассивного участка полета ЛА также используется система (1.1). При этом:

тяга R равна нулю (R=0).

поверхность Земли представляется в виде бесконечной плоскости.

масса ЛА постоянна (m= m-mt=const).

1.3 Алгоритмы численного интегрирования на основе методов Эйлера и Рунге-Кутта

Рассмотрим метод Эйлера — линейное приближение, использующее первые два члена ряда Тейлора. Здесь искомая интегральная кривая аппроксимируется ломаной линией. Если шаг h достаточно мал, то интеграл в формуле можно вычислить используя теорему о среднем, т. е. вынося подынтегральную функцию изпод знака интеграла средним значением.

В методе Эйлера подынтегральная функция выносится при нижнем пределе интегрирования:. Это приближение геометрически соответствует движению от точки x к точке х+h по касательной к кривой y (x) в точке х. Запишем расчетные формулы метода Эйлера:

yk+1=yk+f (xk, yk) h, xk=xk-1+h

y (x0)=y0, yk=y (xk)

В усовершенствованном методе ЭйлераКоши в первом приближении полагается:

а во втором

Погрешность метода Эйлера определяется остаточным членом ряда Тейлора т. е. R~h2 на каждом шаге вычислений. Для обеспечения сходимости шаг h следует выбирать достаточно малым. Для метода ЭйлераКоши погрешность имеет порядок h2.

Рассмотрим метод Рунге и Кутта. В основе получения вычислительных схем этого метода лежит разложение функции y (x) в ряд Тейлора с последующим преобразованием отрезка ряда к виду, не содержащему производных. На шаге h производная dy/dx=f (x, y) аппроксимируется параболой второго порядка. Здесь функция ?(x, h) определяется формулой парабол Симпсона (формула Ньютона — Котеса для трех узлов):

Рассмотрим дифференциальное уравнение при начальном условии (хА, уА). Выполним следующие операции:

По известным начальным условиям (хА, уА) определим значение производной в начальной точке А: .

Из начальной точки, А проведем прямую (рис 1.2.)

и отметим значение ее ординаты в середине шага интегрирования h (точка В с координатами

).

Рис 1.2.

2) Найдем значение производной по формуле в точке В: и проведем из точки, А прямую. Отметим значение ординаты этой прямой в середине шага интегрирования h (точка С с координатами).

3) По уравнению найдем значение производной в точке С: и проведем из точки, А прямую. Отметим значение ординаты этой прямой в конце шага интегрирования h (точка D с координатами).

4) По уравнению найдем значение производной в точке D: .

В результате построений найдем значение производных в точках А, В, С и D. Отложим эти значения на графике рис 1.3. Как видно из графика, в точке с абсциссой получены два значения производной вместо одного. Это следствие приближенности метода. Примем в этой точке среднее значение производной:. Отложив на графике (рис 1.3.) ординату, получим точку М.

Рис 1.3.

Будем считать, что кривая, изображающая зависимость должна проходить через точки A, M и D. Проведем через эти три точки параболу, уравнение которой:

.

Значения коэффициентов a, b и с выбираются из условия прохождения параболы через точки А, М, и D. Коэффициент. Из уравнения параболы имеем систему:

Решив эти уравнения, найдем:

Проинтегрируем теперь уравнение параболы в пределах от x=xA до x=xA+h. Значение этого интеграла является приращением искомой функции y при изменении х на величину h. Таким образом Подставив сюда полученные выше выражения для a, b, c, после приведения подобных членов для общего случая () получим:

Как видно, приращение искомой функции на шаге h при помощи описанных построений удалось представить через значения первых производных функции в четырех точках, лежащих в пределах шага интегрирования h.

Запишем расчетные формулы метода РунгеКутта:

При условии существования у функции производных четвертого порядка погрешность метода является величиной порядка h5.

Для системы дифференциальных уравнений первого порядка данный алгоритм выполняется для каждого уравнения системы параллельно.

1.4 Решение контрольной задачи для одного шага интегрирования методами Эйлера и Рунге-Кутта (Проведение расчетов без ПК)

1.4.1 Начальные условия

Для решения тестового примера на 1 шаг интегрирования системы (1.1) методами Эйлера и РунгеКутта задаем следующие начальные условия:

1.Шаг интегрирования h принимаем равным 0.05 секунды (h=0.05 сек.);

2. V=Vg, где Vgскорость ЛА в момент схода с направляющих, рассчитываемая по формуле:

где

q_g — угол наклона направляющих; q_g =450=0,785 рад;

tgдульное время, рассчитываемые по формуле:

Vg = и tg =, где

=

Vg =

3. q=q_g — угол наклона направляющих; ??g=450=0,785 рад.

4. хначальная абсцисса, рассчитываемая по формуле:

5. уначальная ордината, рассчитываемая по формуле:

1.4.2 Расчет системы уравнений методом Эйлера (1 шаг)

1. Нахождение коэффициентов К1 для каждого уравнения.

По определению коэффициенты К1- это правые части уравнений, вычисленные при начальных условиях.

а. плотность при начальной высоте :

б. лобовое сопротивление при начальной (дульной) скорости:

X=0,5*0.26*0.1 767*1.22 967*(72.9421)2=15.0288

в. Масса при начальном (дульном) времени:

m=60−20*0.1095=57.81

г. К1(1)= (40 000 — 9.81*57.81* 0.785 — 15.0288) / 57.81 = 683.96 м/с2

К1(2)=-9,81*cos450/ 22,708= -0.7761 рад/с К1(3)=22,708*cos450= 3.51 304 м/с К1(4)=22,708*sin450= 3.4857 м/с

2. Результаты одного шага интегрирования по методу Эйлера:

1.4.3 Расчет системы уравнений методом Рунге-Кутта (1 шаг)

1. Нахождение коэффициентов К1 для каждого уравнения.

По определению коэффициенты К1- это правые части уравнений, вычисленные при начальных условиях.

а. плотность при начальной высоте :

б. лобовое сопротивление при начальной (дульной) скорости:

X=0,5*0.26*0.1 767*1.22 967*(72.9421)2=15.0288

в. Масса при начальном (дульном) времени:

m=60−20*0.1095=57.81

г. К1(1)= 8,1 224 487 289

К1(2)= -0,112 446 838

К1(3)= 1,91 383 734

К1(4)= 1,905 149 817

2. Нахождение коэффициентов К2 для каждого уравнения.

По определению коэффициенты К2- это правые части уравнений, вычисленные при полученных новых значениях параметров.

К2(1)= 8,1 419 713 966

К2(2)= -0,9 992 507

К2(3)= 1,241 945 517

К2(4)= 1,2 270 805 529

3. Нахождение коэффициентов К3 для каждого уравнения По определению К3- это правые части уравнений, вычисленные при полученных новых значениях параметров.

К3(1)= 8,1 417 519 898

К3(2)= -0,99 835 329

К3(3)= 1,2 415 239 848

К3(4)= 1,2 282 011 286

4. Нахождение коэффициента К4 для всех уравнений По определению К4- это правые части уравнений, вычисленные в конце шага К4(1)= 8,1 610 849 832

К4(2)= -0,89 854 278

К4(3)= 1,3 930 436 092

К4(4)= 1,364 415 478

5. Результаты одного шага по методу РунгеКутта.

1.4.4 Сравнение решений тестового примера на 1 шаг, полученных вручную (п. 1.4.1 и 1.4.2) и с использованием программного комплекса

Сравнение результатов производится в соответствии с данными, представленными в таблице 1.3.:

Таблица 1.3

Примечание: Результаты компьютерного решения тестового примера взяты из результатов работы программного комплекса.

1.5 Результаты расчета траекторных параметров неуправляемого ЛА

1.5.1 Движение по направляющим

Результатами расчета траекторных параметров данного участка траектории являются значения дульной скорости (vd) и дульного времени (td).

Результаты расчета приведены в таблице 1.4.:

Таблица 1.4

1.5.2 Активный участок траектории ЛА

Расчет траекторных параметров на конце активного участка полета ЛА произведен с помощью программного комплекса с шагом интегрирования равным 0,05

Таблица 1.5

1.5.3 Пассивный участок траектории ЛА

Расчет траекторных параметров на конце пассивного участка полета ЛА произведен с помощью программного комплекса с шагом интегрирования равным 0.05

Таблица 1.6.

1.5.4 Результаты расчета траекторных параметров неуправляемого ЛА для характерных точек траектории

Результаты расчета траекторных параметров неуправляемого ЛА в данном пункте приводятся для характерных точек траектории, показанных на рисунке 1.11.

Рис 1.11.

где Aточка начала движения по направляющим;

Bточка конца движения по направляющим (начало активного участка);

Cточка конца активного участка (начало пассивного участка);

Eнаивысшая точка траектории;

Dточка контакта снаряда с поверхностью;

Результаты расчета представлены в таблице 1.7.

Таблица 1.7

Примечание: Результаты расчета получены при помощи программного комплекса

1.6 Расчет траектории пассивного участка полета ЛА с использованием параболической теории

Параболическая теория изучает движение тяжелой материальной точки в пустоте в стационарном однородном параллельном поле земного тяготения. ЕЕ применяют для ориентировочного расчета траектории ЛА, имеющих малую скорость полета, и для расчета траектории ЛА на больших высотах, где сила сопротивления ничтожна мала.

Основные формулы параболической теории:

Здесь t0=0, y0=yc=228,664 м, x0=xc=270,014 м,

— параметры конца активного участка рассчитанные по данным п 1.5.2

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

Для вершины траектории (точка D) время вычисляется по формуле:

Подставляя tD в формулы 3 и 4, получим X (tD)=6626.495 м Y (tD)=2804.056 м Для точки Е (конца траектории) имеем Y (tE)=0. Тогда Результаты расчетов представлены в сравнении с расчетами реальной траектории в таблице 1.8.

Таблица 1.8

Приложение 1

Графики

1. V (t) методами численного интегрирования и по параболической теории.

2. и (t) методами численного интегрирования и по параболической теории.

3. y (x) (Метод Эйлера и параболическая теория).

4. y (x) (Метод Рунге-Кутта и параболическая теория).

Рис.

Рис.

Рис.

Рис.

Рис.

Рис.

Приложение 2

Текст программы, разработанной для интегрирования методами Эйлера и Рунге-Кутта. (Программа выполнена на Visual Basic 6.0 Professional Edition)

Рис.

Option Base 1

Const r0 = 1.23

Dim r As Double

Dim y As Double

Dim x As Double

Dim V As Double

Dim Q As Double

Dim Xx As Double

Dim m As Double

Dim Yg As Double

Dim Xg As Double

Dim Cx As Double

Dim Vg As Double

Dim S As Double

Dim Mn As Double

Dim Rashod As Double

Dim Tg As Double

Dim Rs As Double

Dim T As Double

Dim h As Double

Dim Qg As Double

Dim i As Double

Dim K (4, 4) As Double

Dim str As String

Dim j As Double

Private Sub Command1_Click ()

lbl.ForeColor = &HFF&

lbl.Caption = «идет расчет»

DoEvents

Yg = CDbl (Text4.Text) * Sin (CDbl (Text3.Text))

Xg = CDbl (Text4.Text) * Cos (CDbl (Text3.Text))

Cx = CDbl (Text5.Text)

Vg = CDbl (Text1.Text)

S = CDbl (Text6.Text)

Mn = CDbl (Text17.Text)

Rashod = CDbl (Text18.Text)

Tg = CDbl (Text2.Text)

Rs = CDbl (Text19.Text)

T = CDbl (Text16.Text)

h = CDbl (Text7.Text)

Qg = CDbl (Text3.Text)

Text33.Text = «»

Text34.Text = «»

Text36.Text = «»

Text40.Text = «»

Text38.Text = «»

'————————————————-;

'Активный участок траектории полета

'Эйлер

x = Xg

y = Yg

V = Vg

Q = Qg

j = Tg

For i = Tg To T Step h

r = r0 * Exp (-y / 7800)

Xx = 0.5 * Cx * r * S * V ^ 2

m = Mn — Rashod * i

V = V + h * (Rs — 9.81 * m * Sin (Q) — Xx) / m

Q = Q — h * 9.81 * Cos (Q) / V

x = x + h * V * Cos (Q)

y = y + h * V * Sin (Q)

If Option1. Value = True And Round (i, 2) = Round (j, 2) Then

j = j + 1

Text33.Text = Text33. Text & Round (i, 2) & vbCrLf

Text34.Text = Text34. Text & Round (V, 1) & vbCrLf

Text36.Text = Text36. Text & Round (Q, 3) & vbCrLf

Text40.Text = Text40. Text & Round (x, 0) & vbCrLf

Text38.Text = Text38. Text & Round (y, 0) & vbCrLf

End If

If i = Tg Then

str = «Значения коэфициентов» & vbCrLf & «на первом шаге Э» & vbCrLf

str = str & «V = «& Round (V, CLng (Text21.Text)) & vbCrLf

str = str & «Q = «& Round (Q, CLng (Text21.Text)) & vbCrLf

str = str & «X = «& Round (x, CLng (Text21.Text)) & vbCrLf

str = str & «Y = «& Round (y, CLng (Text21.Text)) & vbCrLf & vbCrLf

End If

Next i

Text8.Text = Round (V, CLng (Text21.Text))

Text9.Text = Round (Q, CLng (Text21.Text))

Text10.Text = Round (x, CLng (Text21.Text))

Text11.Text = Round (y, CLng (Text21.Text))

'——————————————;

'Рунге-Кутт

x = Xg

y = Yg

V = Vg

Q = Qg

j = Tg

For i = Tg To T Step h

r = r0 * Exp (-y / 7800)

Xx = 0.5 * Cx * r * S * V ^ 2

m = Mn — Rashod * i

K (1, 1) = h * (Rs — 9.81 * m * Sin (Q) — Xx) / m

K (1, 2) = -9.81 * h * Cos (Q) / V

K (1, 3) = h * V * Cos (Q)

K (1, 4) = h * V * Sin (Q)

r = r0 * Exp (-(y + 0.5 * K (1, 4)) / 7800)

Xx = 0.5 * Cx * r * S * (V + 0.5 * K (1, 1)) ^ 2

m = Mn — Rashod * (i + h / 2)

K (2, 1) = h * (Rs — 9.81 * m * Sin (Q + K (1, 2) / 2) — Xx) / m

K (2, 2) = -9.81 * h * Cos (Q + K (1, 2) / 2) / (V + K (1, 1) / 2)

K (2, 3) = h * (V + K (1, 1) / 2) * Cos (Q + K (1, 2) / 2)

K (2, 4) = h * (V + K (1, 1) / 2) * Sin (Q + K (1, 2) / 2)

r = r0 * Exp (-(y + 0.5 * K (2, 4)) / 7800)

Xx = 0.5 * Cx * r * S * (V + 0.5 * K (2, 1)) ^ 2

m = Mn — Rashod * (i + h / 2)

K (3, 1) = h * (Rs — 9.81 * m * Sin (Q + K (2, 2) / 2) — Xx) / m

K (3, 2) = -9.81 * h * Cos (Q + K (2, 2) / 2) / (V + K (2, 1) / 2)

K (3, 3) = h * (V + K (2, 1) / 2) * Cos (Q + K (2, 2) / 2)

K (3, 4) = h * (V + K (2, 1) / 2) * Sin (Q + K (2, 2) / 2)

r = r0 * Exp (-(y + K (3, 4)) / 7800)

Xx = 0.5 * Cx * r * S * (V + K (3, 1)) ^ 2

m = Mn — Rashod * (i + h)

K (4, 1) = h * (Rs — 9.81 * m * Sin (Q + K (3, 2)) — Xx) / m

K (4, 2) = -9.81 * h * Cos (Q + K (3, 2)) / (V + K (3, 1))

K (4, 3) = h * (V + K (3, 1)) * Cos (Q + K (3, 2))

K (4, 4) = h * (V + K (3, 1)) * Sin (Q + K (3, 2))

V = V + (K (1, 1) + K (4, 1) + 2 * (K (2, 1) + K (3, 1))) / 6

Q = Q + (K (1, 2) + K (4, 2) + 2 * (K (2, 2) + K (3, 2))) / 6

x = x + (K (1, 3) + K (4, 3) + 2 * (K (2, 3) + K (3, 3))) / 6

y = y + (K (1, 4) + K (4, 4) + 2 * (K (2, 4) + K (3, 4))) / 6

If Option2. Value = True And Round (i, 2) = Round (j, 2) Then

j = j + 1

Text33.Text = Text33. Text & Round (i, 2) & vbCrLf

Text34.Text = Text34. Text & Round (V, 1) & vbCrLf

Text36.Text = Text36. Text & Round (Q, 3) & vbCrLf

Text40.Text = Text40. Text & Round (x, 0) & vbCrLf

Text38.Text = Text38. Text & Round (y, 0) & vbCrLf

End If

If i = Tg Then

str = str & «Значения коэффициентов» & vbCrLf & «на первом шаге Р-К» & vbCrLf

str = str & «K (1,1) = «& Round (K (1, 1), CLng (Text21.Text)) & vbCrLf

str = str & «K (1,2) = «& Round (K (1, 2), CLng (Text21.Text)) & vbCrLf

str = str & «K (1,3) = «& Round (K (1, 3), CLng (Text21.Text)) & vbCrLf

str = str & «K (1,4) = «& Round (K (1, 4), CLng (Text21.Text)) & vbCrLf

str = str & «K (2,1) = «& Round (K (2, 1), CLng (Text21.Text)) & vbCrLf

str = str & «K (2,2) = «& Round (K (2, 2), CLng (Text21.Text)) & vbCrLf

str = str & «K (2,3) = «& Round (K (2, 3), CLng (Text21.Text)) & vbCrLf

str = str & «K (2,4) = «& Round (K (2, 4), CLng (Text21.Text)) & vbCrLf

str = str & «K (3,1) = «& Round (K (3, 1), CLng (Text21.Text)) & vbCrLf

str = str & «K (3,2) = «& Round (K (3, 2), CLng (Text21.Text)) & vbCrLf

str = str & «K (3,3) = «& Round (K (3, 3), CLng (Text21.Text)) & vbCrLf

str = str & «K (3,4) = «& Round (K (3, 4), CLng (Text21.Text)) & vbCrLf

str = str & «K (4,1) = «& Round (K (4, 1), CLng (Text21.Text)) & vbCrLf

str = str & «K (4,2) = «& Round (K (4, 2), CLng (Text21.Text)) & vbCrLf

str = str & «K (4,3) = «& Round (K (4, 3), CLng (Text21.Text)) & vbCrLf

str = str & «K (4,4) = «& Round (K (4, 4), CLng (Text21.Text)) & vbCrLf

str = str & «V = «& Round (V, CLng (Text21.Text)) & vbCrLf

str = str & «Q = «& Round (Q, CLng (Text21.Text)) & vbCrLf

str = str & «X = «& Round (x, CLng (Text21.Text)) & vbCrLf

str = str & «Y = «& Round (y, CLng (Text21.Text))

Text20.Text = str

End If

Next i

Text12.Text = Round (V, CLng (Text21.Text))

Text13.Text = Round (Q, CLng (Text21.Text))

Text14.Text = Round (x, CLng (Text21.Text))

Text15.Text = Round (y, CLng (Text21.Text))

'————————————————-;

'Пассивный участок траектории полета

'Эйлер

V = CDbl (Text8.Text)

Q = CDbl (Text9.Text)

x = CDbl (Text10.Text)

y = CDbl (Text11.Text)

i = T

m = Mn — Rashod * T

j = T

Do While y > 0

r = r0 * Exp (-y / 7800)

Xx = 0.5 * Cx * r * S * V ^ 2

V = V — h * (9.81 * m * Sin (Q) + Xx) / m

Q = Q — h * 9.81 * Cos (Q) / V

x = x + h * V * Cos (Q)

y = y + h * V * Sin (Q)

If Option1. Value = True And Round (i, 2) = Round (j, 2) Then

j = j + 1

Text33.Text = Text33. Text & Round (i, 2) & vbCrLf

Text34.Text = Text34. Text & Round (V, 1) & vbCrLf

Text36.Text = Text36. Text & Round (Q, 3) & vbCrLf

Text40.Text = Text40. Text & Round (x, 0) & vbCrLf

Text38.Text = Text38. Text & Round (y, 0) & vbCrLf

End If

If Abs (Q) < 0.001 Then

str = «Значения параметров» & vbCrLf & «в наивысшей точке траектории» & vbCrLf & «методом Эйлера» & vbCrLf

str = str & «t = «& Round (i, 2) & vbCrLf

str = str & «V = «& Round (V, CLng (Text21.Text)) & vbCrLf

str = str & «Q = 0» & vbCrLf

str = str & «X = «& Round (x, CLng (Text21.Text)) & vbCrLf

str = str & «Y = «& Round (y, CLng (Text21.Text)) & vbCrLf & vbCrLf

End If

i = i + h

Loop

Text25.Text = Round (V, CLng (Text21.Text))

Text24.Text = Round (Q, CLng (Text21.Text))

Text23.Text = Round (x, CLng (Text21.Text))

Text22.Text = Round (y, CLng (Text21.Text))

Text31.Text = Round (i, 2)

'——————————————;

'Рунге-Кутт

V = CDbl (Text12.Text)

Q = CDbl (Text13.Text)

x = CDbl (Text14.Text)

y = CDbl (Text15.Text)

i = T

j = T

m = Mn — Rashod * T

Do While y > 0

r = r0 * Exp (-y / 7800)

Xx = 0.5 * Cx * r * S * V ^ 2

K (1, 1) = h * (-9.81 * m * Sin (Q) — Xx) / m

K (1, 2) = -9.81 * h * Cos (Q) / V

K (1, 3) = h * V * Cos (Q)

K (1, 4) = h * V * Sin (Q)

r = r0 * Exp (-(y + 0.5 * K (1, 4)) / 7800)

Xx = 0.5 * Cx * r * S * (V + 0.5 * K (1, 1)) ^ 2

K (2, 1) = h * (-9.81 * m * Sin (Q + K (1, 2) / 2) — Xx) / m

K (2, 2) = -9.81 * h * Cos (Q + K (1, 2) / 2) / (V + K (1, 1) / 2)

K (2, 3) = h * (V + K (1, 1) / 2) * Cos (Q + K (1, 2) / 2)

K (2, 4) = h * (V + K (1, 1) / 2) * Sin (Q + K (1, 2) / 2)

r = r0 * Exp (-(y + 0.5 * K (2, 4)) / 7800)

Xx = 0.5 * Cx * r * S * (V + 0.5 * K (2, 1)) ^ 2

K (3, 1) = h * (-9.81 * m * Sin (Q + K (2, 2) / 2) — Xx) / m

K (3, 2) = -9.81 * h * Cos (Q + K (2, 2) / 2) / (V + K (2, 1) / 2)

K (3, 3) = h * (V + K (2, 1) / 2) * Cos (Q + K (2, 2) / 2)

K (3, 4) = h * (V + K (2, 1) / 2) * Sin (Q + K (2, 2) / 2)

r = r0 * Exp (-(y + K (3, 4)) / 7800)

Xx = 0.5 * Cx * r * S * (V + K (3, 1)) ^ 2

K (4, 1) = h * (-9.81 * m * Sin (Q + K (3, 2)) — Xx) / m

K (4, 2) = -9.81 * h * Cos (Q + K (3, 2)) / (V + K (3, 1))

K (4, 3) = h * (V + K (3, 1)) * Cos (Q + K (3, 2))

K (4, 4) = h * (V + K (3, 1)) * Sin (Q + K (3, 2))

V = V + (K (1, 1) + K (4, 1) + 2 * (K (2, 1) + K (3, 1))) / 6

Q = Q + (K (1, 2) + K (4, 2) + 2 * (K (2, 2) + K (3, 2))) / 6

x = x + (K (1, 3) + K (4, 3) + 2 * (K (2, 3) + K (3, 3))) / 6

y = y + (K (1, 4) + K (4, 4) + 2 * (K (2, 4) + K (3, 4))) / 6

If Option2. Value = True And Round (i, 2) = Round (j, 2) Then

j = j + 1

Text33.Text = Text33. Text & Round (i, 2) & vbCrLf

Text34.Text = Text34. Text & Round (V, 1) & vbCrLf

Text36.Text = Text36. Text & Round (Q, 3) & vbCrLf

Text40.Text = Text40. Text & Round (x, 0) & vbCrLf

Text38.Text = Text38. Text & Round (y, 0) & vbCrLf

End If

If Abs (Q) < 0.001 Then

str = str & «Значения параметров» & vbCrLf & «в наивысшей точке траектории» & vbCrLf & «методом Рунге-Кутта» & vbCrLf

str = str & «t = «& Round (i, 2) & vbCrLf

str = str & «V = «& Round (V, CLng (Text21.Text)) & vbCrLf

str = str & «Q = 0» & vbCrLf

str = str & «X = «& Round (x, CLng (Text21.Text)) & vbCrLf

str = str & «Y = «& Round (y, CLng (Text21.Text))

Text30.Text = str

End If

i = i + h

Loop

Text29.Text = Round (V, CLng (Text21.Text))

Text28.Text = Round (Q, CLng (Text21.Text))

Text27.Text = Round (x, CLng (Text21.Text))

Text26.Text = Round (y, CLng (Text21.Text))

Text32.Text = Round (i, 2)

lbl.ForeColor = &HC000&

lbl.Caption = «Расчет выполнен»

End Sub

Private Sub Command2_Click ()

End

End Sub

Список литературы

1. Бахвалов Н. С. Численные методы. т.1., М.: Наука, 1975.

2. Волков Е. А. Численные методы. М.: Наука, 1982

3. Шалыгин А. С. и др. Автоматизация расчета траекторий ЛА: Учебн. пособие. — Л., 1990.

4. http://a5.bstu.spb.su

Размещено на