Дифференциальных уравнений
Ряд задач, возникающих перед инженерами-механиками при проектировании или исследованиях металлургических машин, сводится к решению обыкновенных дифференциальных уравнений. С помощью таких уравнений описываются модели динамики объектов, как правило, они отражают изменение параметров объектов во времени.
В общем случае обыкновенное дифференциальное уравнение имеет вид:
или y¢ = f(x,y) . (65)
Решением дифференциального уравнения является функция y= y(x), отыскание которой выполняется интегрированием дифференциального уравнения. При использовании численных методов решение уравнения (65) представляется в табличном виде, т.е. получается совокупность значений yi и xi.
Одним из основных типов задач для обыкновенных дифференциальных уравнений является задача Коши, которая для дифференциального уравнения первого порядка заключается в нахождении решения уравнения:
y¢ = f (x,y) , (66)
удовлетворяющего при x=x0 начальному условию
y(x0) = y0 . (67)
Интегрируя (66) в пределах от x до x+h (где h – шаг), получим:
. (68)
В зависимости от вида функции, аппроксимирующей подынтегральную функцию (константа; прямая, проходящая через концы отрезка; парабола и т.д.) получают различные формулы численного интегрирования дифференциального уравнения.
Решение дифференциального уравнения численными методами носит шаговый характер, т.е. по одной или по нескольким начальным точкам (x, y) за один шаг находят следующую точку, затем следующую и т.д. Разница между двумя соседними значениями аргумента x называется шагом: h = xi+1 - xi.
Выделяют два класса методов решения: одношаговые и многошаговые. Первый класс методов требует для нахождения следующего значения функции y= y(x) только одной текущей точки, а второй – нескольких. Поэтому методами второго класса нельзя начать решение задачи Коши, это всегда делается одношаговыми методами.
Основная идея получения простейших вычислительных алгоритмов в одношаговых методах сводится к следующему. Зная значение y в точке разложения yi и производную f(xi , yi), находят значение функции y(х) через шаг h:
yi+1 = yi + Dyi. (69)
Основными одношаговыми методами решения дифференциальных уравнений являются методы Эйлера и Рунге-Кутта [6, 7]. Рассмотрим их вычислительные алгоритмы.
Использование любого из методов начинается с выбора фиксированного шага h по переменной х:
, (70)
где xк – конечная точка интервала интегрирования;
n - количество шагов интегрирования.
Далее используется одна из приближенных формул вычисления функции y=y(x) в точке x0+h. Принимая эту точку в качестве начальной, определяют y(x0+2h). Этот прием повторяют необходимое количество раз и вычисляют значения y(x) в конечной точке. При необходимости сохраняют массив значений y(x0+ih), что дает возможность построить график y=y(x).
Наиболее простым методом решения дифференциальных уравнений первого порядка (66) является метод Эйлера.
Для определения значений функции y=y(x) используется выражение:
yi+1 = yi + h f(xi, yi) , i = 0, 1, .... . (71)
Для начала вычислений используются начальные условия (67).
На рисунке 35 приведена графическая интерпретация метода Эйлера. В точке (xi , yi) функция y=y(x) аппроксимируется отрезком, проходящим через нее под углом ai. Угол ai определяется значением производной в точке (xi , yi), т.е. tg ai = f(xi , yi).
Рисунок 35 – Графическая интерпретация метода Эйлера
Блок-схема алгоритма решения представлена на рисунке 36.
Рисунок 36 – Блок-схема алгоритма решения дифференциального уравнения методом Эйлера
Метод Эйлера обладает большой погрешностью и имеет систематическое накопление ошибок (см. рисунок 35). Погрешность метода пропорциональна h2.
Для повышения точности на практике используют модифицированный метод Эйлера. Он имеет следующий вычислительный алгоритм:
yi+1 = yi + h f(xi + 0.5h , yi + 0.5h f(xi, yi)) , i = 0, 1, .... . (72)
На рисунке 37 изображено графическое пояснение формулы (72). Здесь через точку (xi , yi) проводится прямая L1 под углом ai ( tg ai = f(xi , yi) ) и определяется точка А с координатами (xi+0.5h, yi+0.5h f(xi, yi)). В точке А вычисляется производная, которая определяет угол наклона прямой L2 , т.е. tg bi = f(xi+0.5h, yi + 0.5h f(xi, yi)). Через точку (xi , yi) параллельно L2 проводится прямая L до пересечения с вертикальной прямой, исходящей из точки на оси х с абсциссой xi+1 . Отсекаемый при этом отрезок на вертикальной прямой и определяет значение функции y=y(x) равное yi+1 при x=xi+1, соответствующее вычисленному по (72).
Рисунок 37 – Графическое пояснение вычислений по модифицированному методу Эйлера
Блок-схема алгоритма решения дифференциального уравнения модифицированным методом Эйлера приведена на рисунке 38. Этот метод дает погрешность пропорциональную h3.
Наиболее распространенным методом решения задачи Коши (66)-(67) является метод Рунге-Кутта четвертого порядка. Этот метод более точен по сравнению с методами Эйлера. Для расчета одного значения функции необходимо четыре раза вычислить правую часть дифференциального уравнения (66), а не два, как в модифицированном методе Эйлера.
Рисунок 38 – Блок-схема алгоритма решения дифференциального уравнения модифицированным методом Эйлера
Вычислительный алгоритм метода основывается на следующих выражениях:
,
k1 =h f(xi , yi) , k2 =h f(xi+0.5h , yi+0.5k1) , (73)
k3 =h f(xi+0.5h , yi+0.5k2) , k4 =h f(xi+h , yi+ k3) , i = 0, 1, .... .
Блок-схема алгоритма решения дифференциального уравнения методом Рунге-Кутта четвертого порядка приведена на рисунке 39.
Погрешность метода Рунге-Кутта пропорциональна h5.
Рисунок 39 – Блок-схема алгоритма решения дифференциального уравнения методом Рунге-Кутта четвертого порядка
СПИСОК ЛИТЕРАТУРЫ
1. Тимошенко Г.М., Зима П.Ф. Теория инженерного эксперимента. – К.: УМК ВО, 1991. – 124 с.
2. Боглаев Ю.П. Вычислительная математика и программирование. – М.: Высш. шк., 1990. – 544 с.
3. Дьяконов В.П. Справочник по алгоритмам и программам на языке Бейсик для персональных ЭВМ. – М.: Наука, 1987. – 240 с.
4. Мудров А.Е. Численные методы для ПЭВМ на языках Бейсик, Фортран и Паскаль. – Томск: Раско, 1991. – 272 с.
5. TURBO PASCAL: Алгорітми і програми: Чисельні методи в фізиці та математиці. – К.: Вища шк., 1992. – 247 с.
6. Васильков Ю.В., Василькова Н.Н. Компьютерные технологии вычислений в математическом моделировании. – М.: Финансы и статистика, 1999. – 256 с.
7. Гловацкая А.П. Методы и алгоритмы вычислительной математики. – М.: Радио и связь, 1999. – 408 с.
8. Пироженко Н.Г. Расчет параметров привода металлургических машин. – Донецк: ДонГТУ: кафедра МОЗЧМ. – 1999. – 84 с.
9. Левин М.З., Седуш В.Я. Механическое оборудование доменных цехов. – Киев-Донецк: Вища школа, 1978. – 1978. – 176 с.
10. Механическое оборудование сталеплавильных цехов / Левин М.З., Седуш В.Я., Мачикин В.И. и др. – Киев-Донецк: Вища школа, 1985. – 165 с.