НОУ ИНТУИТ | Введение в математическое моделирование. Лекция 12: Компьютерное моделирование и решение нелинейных уравнений

Учитесь и получайте официальные документы БЕСПЛАТНО. Вы можете поддержать наш проект.

Регистрация Вход

Твой путь к знаниям!

Опубликован: 15.03.2007 | Уровень: специалист | Доступ: платный | ВУЗ: Донецкий национальный технический университет

|

Вам нравится? Нравится 69 студентам

| Поделиться |

Поддержать программу

Решение дифференциальных уравнений m-го порядка методом Рунге-Кутта (4-го порядка)

Как уже было сказано, любое дифференциальное уравнение m-го порядка

$y_{(m)}=f(x,y,y',y'',\ldots,y^{(m-1)})$

( 12.15)

сводится к системе, состоящей из m дифференциальных уравнений 1-го порядка

$\left\{ \begin{array}{l} y' = y_1\\ y'_1 = y_2\\ y'_2 = y_3\\ \ldots\\ y'_{(m-1)} = f(x,y_1,y_2,\ldots,y_{(m-1)} \end{array} \right.$

( 12.6)

Численным решением системы (12.9), а значит и дифференциального уравнения m-го порядка (12.8) является m табличных функций

$y, y_1=y', y_2=y'', \ldots, y_m=y^{(m-1)},$

т.е. функция y(x) и все ее производные, включая производную (m-1)-го порядка.

При этом каждая из табличных функций определяется на промежутке [a, b] с шагом h и включает n узловых точек. Таким образом, численным решением уравнения (12.8) или системы (12.9) является матрица порядка $(n\times m)$ , (табл. 12.8)

Таблица 12.1.
i	x	y	y₁=y',	y₁=y''₁,	:	y_m-1=y^(m-1)
0	x₀	y₀	(y₁)₀	(y₂)₀	:	(y_m-1)₀
1	x₁	y₁	(y₁)₁	(y₂)₁	:	(y_m-1)₁
2	x₂	y₂	(y₁)₂	(y₂)₂	:	(y_m-1)₂
3	x₃	y₃	(y₁)₃	(y₂)₃	:	(y_m-1)₃
:	:	:	:	:	:	:
n	x_n	y_n	(y₁)_n	(y₂)_n	:	(y_m-1)_n

где

m - порядок дифференциального уравнения, равен количеству столбцов матрицы,

n = (b-a)/h - количество шагов интегрирования, равно количеству строк матрицы.

Каждый j -й столбец матрицы - это массив решений одной j -й табличной функции по всем n шагам интегрирования.

Каждая i -ая строка матрицы - это массив решений m табличных функций на одном i -ом шаге интегрирования.

На графике решением дифференциального уравнения m -го порядка (12.8) является совокупность $(n\times m)$ узловых точек. При этом каждому шагу интегрирования, т.е. каждому значению x_i, $i=\overline{1,n}$ , соответствуют m узловых точек скоординатами

$(x_i,y_1), (x_i,(y_1)_i), (x_i,(y_2)_i), \ldots, (x_i,(y_{m-1})_i).$

При построении алгоритма задачи будем как и ранее расчет вести по шагам интегрирования, т.е. в цикле по $i=\overline{1,n}.$ При этом, как и ранее, на каждом i -ом шаге цикла будем рассчитывать решение дифференциального уравнения и тут же его печатать. Тогда нет необходимости формировать матрицу решений, а можно ограничиться формированием массива решений (12.15), который соответствует одной i-й строке матрицы:

$Y=[y(1), y(2), y(3), \ldots , y(m)],$

( 12.17)

где

$y(1)= y(x),\\ y(2)= y'(x),\\ y(3)= y"(x),\\ \ldots \ldots\\ y(m)= y^{m-1}(x).$

( 12.18)

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

Рассмотрим пример (12.19).

Дано дифференциальное уравнение второго порядка

$y''=\frac{1}{x}y' - \frac{x^2-1}{x^2}y$

( 12.19)

С учетом обозначений (12.18) имеем:

$y(1)=y(x),\\ y(2)=y'(x)$

Тогда дифференциальное уравнение (12.19) сводится к системе:

$\left\{ \begin{array}{l} y'(1)=y(2),\\ y'(2)=\frac{1}{x}y(2)-\frac{x^2-1}{x^2}y(1). \end{array} \right.$

( 12.20)

Вычисление правых частей уравнений системы (12.20) будем выполнять в подпрограмме PRAV. При этом подпрограмма PRAV будет иметь вид

$\left\{ \begin{array}{l} f(1)=y(2),\\ f(2)=\frac{1}{x}y(2)-\frac{x^2-1}{x^2}y(1). \end{array} \right.$

Обращение к подпрограмме

PRAV(m, x, Y, F),

где

m -порядок системы,

x - значение x на i -м шаге интегрирования,

Y=[y(1), y(2), y(3), : , y(m)] - входной массив длиной m,

F=[f(1), f(2),:, f(m)] - результат, массив значений правых частей уравнений системы (12.10) длиной m.

Программу решения системы дифференциальных уравнений реализуем в виде 3-х программных модулей:

основной программы, в которой организуем циклический процесс по всем шагам интегрирования, $i=\overline{1,n}$ ;
подпрограммы RGK, в которой на каждом i -м шаге реализуется метод Рунге Кутта (4-го порядка) для системы дифференциальных уравнений m -го порядка. Здесь на каждом шаге интегрирования вычисляется массив решений Y длиной m. Для этого внутри подпрограммы RGK организуем циклы по j, где $j=\overline{1,m}$ ;
подпрограммы PRAV, обращение к которой осуществляется из подпрограммы RGK для вычисления массива F длиной m - значений правых частей уравнений системы.

Схема алгоритма основной программы представлена на рис.12.17.

Рис. 12.17. Схема алгоритма основной программы

Здесь

m -порядок системы,

h -шаг интегрирования,

n -количество шагов интегрирования,

x -начальное и далее - текущее значение x,

Y -массив длиной m, куда заносим начальные и далее - текущие значения решений системы на одном шаге интегрирования.

В подпрограмме RGK для вычисления элементов массива Y, используем те же формулы, что и для решения одного дифференциального уравнения 1-го порядка методом Рунге-Кутта (4-го порядка), но с учетом поправки на массивы.

Тогда

$T(1,j)=h\cdot f_ j(x, Y),\\ y1(j)=y(j)+T(1,j)/2, j=\overline{1,m},\\ T(2,j)=h\cdot f_ j(x+h/2, Y1),\\ y1(j)=y(j)+T(2,j)/2, j=\overline{1,m},\\ T(3,j)=h\cdot f_ j(x+h/2, Y1),\\ y1(j)=y(j)+T(3,j)/2, j=\overline{1,m},\\ T(4,j)=h\cdot f_ j(x+h, Y1),\\ y(j)=y(j)+\frac{1}{6}(T(1,j)+2T(2,j)+2T(3,j)+T(4,j)).$

Здесь

Y - массив решений длиной m,

Y1 - рабочий массив длиной m,

T - рабочая матрица порядка ( $4\times m$ ).

Для вычисления значений f_j правых частей дифференциальных уравнений системы перед вычислением каждого элемента матрицы Т на каждом j -ом шаге необходимо выполнить обращение к подпрограмме PRAV.

Схема алгоритма подпрограммы RGK представлена на рис.12.18

Рис. 12.18. Схема алгоритма подпрограммы RGK

Дальше >>

Введение в математическое моделирование

Введение в математическое моделирование

Компьютерное моделирование и решение нелинейных уравнений

Решение дифференциальных уравнений m-го порядка методом Рунге-Кутта (4-го порядка)

Вопросы и ответы

Студенты

Авторизоваться

Введение в математическое моделирование

Введение в математическое моделирование

Компьютерное моделирование и решение нелинейных уравнений

Решение дифференциальных уравнений m-го порядка методом Рунге-Кутта (4-го порядка)

Вопросы и ответы

Студенты