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

 

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

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

y_k+1=y_k+f(x_k,y_k)h, x_k=x_k-1+h

y(x₀)=y₀, y_k=y(x_k)

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

а во втором

Погрешность метода Эйлера определяется остаточным членом ряда Тейлора

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

 

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

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

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

Из начальной точки А проведем прямую (рис 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=x_A до x=x_A+h. Значение этого интеграла является приращением искомой функции y при изменении х на величину h. Таким образом Подставив сюда полученные выше выражения для a,b,c, после приведения подобных членов для общего случая () получим:

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

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

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

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

 

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

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

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

 

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

 

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

, где

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

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

V_g = и t_g = , где

=

V_g =

3. q=q_g- угол наклона направляющих; q_g=45⁰=0,785 рад.

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

х=S_n*cosq_g=4*0.785=3.96 м

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

у= S_n*sinq_g=4*0.785=3.96 м

 

 

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