Методы анализа нелинейных динамических моделей. - Холодниок М.
Скачать (прямая ссылка):
dz
где коэффициенты ^ подсчитываются методом исключения Гаусса. Выбор
индекса k представляет собой самостоятельную проблему. Мы либо фиксируем
его заранее, либо видоизменяем алгоритм исключения Гаусса с выбором
главного элемента для решения (неквадратной) системы п линейных
алгебраических уравнений (5.2.12) относительно п+1 неизвестных. При этом
неизвестная, в столбце которой отсутствует главный элемент, есть
неизвестная хн
зависимости
Подставим (5.2.9);
найденные
(5.2.13) в уравнение
dx
k
dz
= ±
[я + 1 1 !/2
¦¦и
i?>k J
(5.2.14)
*> Речь идет о выборе главного элемента "по строке". При этом переменные
могут быть предварительно масштабированы, см. [5.10]. - Прим. ред.
5.2. Зависимость стационарных решений от параметра
141
Таким образом, мы определили производную dxk/dz с точностью до знака.
Если индекс k остается неизменным для всего продолжения, то знак в
формуле (5.2.14) в процессе продолжения не меняется *>. Отсюда следует,
что тогда хк во время всего процесса продолжения должно либо постоянно
возрастать, либо постоянно убывать. Однако это не всегда бывает так, и
тогда фиксированный выбор k невозможен. Поэтому в практических ситуациях
выбор k осуществляется адаптивным образом-- значение индекса k получается
в результате применения метода исключения Гаусса (см. выше). При этом
знак в формуле (5.2.14) будет определяться тем, возрастает или убывает в
этот момент выбранная величина хк вдоль кривой в направлении возрастания
z. Поэтому мы введем "знаковые" переменные Ni, положив
Ni = sign(^-^, t = 1, ..., я + 1. (5.2.15)
Указанные "параметры направления" будут вычисляться на каждом шаге
процесса продолжения. Выбор знака в формуле
(5.2.14) определяется, следовательно, соответствующей величиной Nk,
взятой с предыдущего шага. Отметим, что оставшиеся производные dxi/dz
(i^k) вычисляются по формуле (5.2.13), где коэффициенты [3; известны.
Затем после вычисления этих производных мы вновь вычисляем параметры
направления Ni и делаем шаг по параметру г. Соответствующая расчетная
схема представлена на рис. 5.3, причем для интегрирования
дифференциальных уравнений (5.2.13) и (5.2.14) в данном случае
используется обычный метод Эйлера.
Проиллюстрируем использование описанного метода продолжения на примере
нахождения зависимости стационарного решения задачи 1 от времени задержки
т (ср. рис. 5.2). Уравнения (5.1.3) в указанном случае переписываются в
виде
/, (х, 0, т) = -Ах + V (1 - х) ехр = О,
f2(x, 0, т) = - Л0 + k0xB(l - .t)expTq^- - ат(0 - 0С) = О.
(5.2.16)
Начальное условие (5.2.10) при х - 0 принимает вид z = 0: х = 0, 0 =
0, т = 0,
dxk
11 Пока J* невырожденна, |3; конечны и - не обращается в 0.- Прим. ред.
142
Глава 5
Рис. 5.3. Схема метода дифференцирования по длине дуги и метода
продолжения Эйлера.
Таблица 5.4. Результаты вычислений по методу продолжения с помощью
алгоритма, представленного на рис. 5.3 (задача 1: у = 20, В = 10, Л = I,
0 = -5, k0 = 1, а=1). Исходные значения коэффициентов направления: Ni =
1, Л12 = 1, Nз = 1, масштабные коэффициенты в методе Гаусса 4)
Р1 = 1> ?2 - 0,1, Ръ = 1, Аг = 0,01.
Z X е •с II 111 N, N,
0 0,0000 0,0000 0,0000 0,0 1 1 1
0,5 0,0780 0,4905 0,0527 4,9 Е - 3 1 1 1
1 0,1355 0,9870 0,0615 7,5Е - 3 1 1 1
1,5 0,1857 1,484 0,0575 8,6Е - 3 1 1 1
2,0 0,2333 1,982 0,0503 8,9Е - 3 1 1 J
3,0 0,3268 2,978 0,0364 7,9Е - 3 1 1 1
4,0 0,4210 3,973 0,0265 5,7Е - 3 1 1 - 1
5,0 0,5168 4,969 0,0200 2,9 Е - 3 1 1 - 1
7,0 0,7129 6,959 0,0142 4,6Е - 3 1 1 - 1
8,0 0,8145 7,954 0,0148 1,ЗЕ - 2 1 1 1
9,0 0,9366 8,946 0,0302 1,2Е - 1 1 1 1
9,5 0,9743 8,516 0,0908 5,8Е - 1 1 - 1 1
10,0 0,9778 8,018 0,1352 6,0Е - 1 1 - 1 1
11,0 0,9779 7,022 0,2292 6,1Е - 1 - 1 - 1 1
13,0 0,9649 5,036 9,4592 6,ЗЕ - 1 - 1 - 1 1
16,0 0,8690 2,076 0,9336 6,7Е - 1 - 1 - 1 1
18,0 0,5883 0,1114 1,127 6,9 Е - 1 - 1 - 1 - 1
19,0 0,3483 -0,8555 1,044 6.9 Е - 1 - 1 - 1 - 1
20,0 0,1980 -1,825 1,195 6,9 Е - 1 - 1 - 1 1
21,0 0,1381 -2,663 1,724 6,8Е - 1 - 1 - 1 1
Д z
0,2 0,1 0,05 '0,02 0,01 0,005
г 5 5 5 5 5 5
X 0,5163 0,5166 0,5167 0,5168 0,5168 0,5168
0 4,9664 4,9676 4,9682 4,9685 4,9686 4,9686
т 0,0201 0,0201 0,0200 0,0200 0,0200 0,0200
Ilf II 2,9 Е - 2 2,2Е - 2 1,ЗЕ - 2 5,6Е - 3 2,9Е - 3 1.5Е-
3
2 12 12 12 12 12 12
X 0,9665 1,1904 0,9972 0,9997 0,9736 0,9797
0 11,814 11,912 11,869 11,818 6,0281 6,3726
т -0,1286 -0,0005 -0,1128 -0,1084 0,3360 0,3010
Ilf II 8,2Е1 1.0Е1 1.5Е1 1.1Е1 6,2Е - 1 2ДЕ0
1; При выборе "главного" элемента в строке 1 берется max (р; а,у).
Аналогично на следующих шагах метода исключения.
144
Глава 5
а (расширенная) матрица Якоби J системы (5.2.16) (если обозначить Е -
ехр(0/(1 + 0/v))) записывается как
Применим теперь вычислительный алгоритм, соответствующий расчетной схеме
рис. 5.3, для случая различных шагов интегрирования Да. Результаты
вычислений представлены в табл. 5.4. Первая часть этой таблицы может быть
использована для более глубокого уяснения метода. Читатель может легко
