Слаботочка Книги

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 [132] 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208

В частном случае А,Ъ = const имеем u(i) = u){t)b, где матрица uj{t) имеет вид

= exp{A{t - т)} dr = (ехр {At} - .

Значение ш{Н) можно было бы попытаться вычислить, пользуясь разложением в ряд Тейлора

uj{H)HY-{AHy~\ (9)

Однако в случае жестких систем при реально допустимых значениях Н величина ИЯ 1 и

а) для достижения приемлемой точности требуется взять слишком много слагаемых;

б) в тех случаях, когда требуемое для достижения нужной точности число слагаемых допустимо по числу арифметических действий, использование разложения (9) при ЦЛЦЯ > 1 может быть неприемлемо по другой причине: среди учитываемых в правой части слагаемых встречаются очень большие, и относительная погрешность, вызванная округлениями, недопустимо велика.

Матрица a)(i) удовлетворяет соотношению

-W=-()(2i + Q)), (10)

где Е - единичная матрица. Поэтому часто бывает целесообразно пойти по следующему пути. Выбираем s такое, что ЛЯ-2 * < 1; основываясь на (9), вычисляем

V2V ~ 2 г! V 2

а затем ш{Н/2 ),..., u)(H/2), и){Н) с помощью рекуррентной формулы (10).

В случае линейной системы у = А{х)у алгоритм решения задачи может быть несколько упрощен. Здесь в описанном выше алгоритме на каждом шаге возникает необходимость решения системы у = А{хп)у при начальном условии у(ж ) = у . Тогда

у +1 = ip{H)yn, где ip{H) = ехр {АН}, А = Л(ж ).

На первый взгляд может показаться разумным следующий путь. Матрица ip{H) удовлетворяет соотношению

а затем

пользуясь рекуррентной формулой (И). Такой путь при большом л-приводит к суш,ественному накоплению погрешности. Поэтому положим ф{Н) = <р{Н) - Е, вычислим

и затем ?/(i 2*~),..., г/) (iJ), пользуясь рекуррентной формулой Далее находим ip{H.) - Е Л- г\){Щ.

ф{г) = ф(~] 2Е + ф{

При {х ф о метод точности 0{fi) получается, если в (7) также учесть второе слагаемое. Тогда потребуется аналогичным образом сконструировать точный метод решения вспомогательного уравнения

U = ли + b + cf, с = --

В случае решения линейной задачи и = А{х)и + f(x) можно предложить довольно простые методы точности О (/).*).

2. Другая группа методов решения жестких задач строится следуюпщм образом. Зададимся некоторым к и приблизим производную у(ж ) односторонней аппроксимацией к-го порядка точности

-г=1 - г=0

Выражение Г(ж , у(ж )) оставим без изменения. Получим конечно-разностную аппроксимацию

E-f(n,yn)=0. (12)

г=0

Рассмотрим случай модельного уравнения у = My, когда (12) превращается в конечно-разностное уравнение к

а-гУп-i

Е-му = о.

поэтому зададимся некоторым s и вычисляем

1 - M/i- 1 - Mh Если Re М < О, в частности, если М действительно и М < О, то при любом h < 1 и погрешности, возникающие в процессе счета, убывают.

l~Mh

Поэтому такой метод может быть применен к решению жестких систем.

Решение этого уравнения вьшисывается через корни характеристического уравнения

a-i/i- - /гМ/х = 0.

При /с = 1, 2 корни этого уравнения удовлетворяют условию 1 в области значений М : ReM < 0; при к = 3, 4, 5, 6- условию 1 в области значений М : 11ш М -а. Re М, где Ок > 0.

В нелинейном случае значение у требуется находить из нелинейной системы (12).

А.пгоритмы решения жестких систем различаются способами нахождения начального приб.пижения к решению (12) и алгоритмами приближенного решения (12). Рассмотрим простейший случай к - 1, когда (12) превращается в неявный метод Эйлера:

Могло бы показаться разумным найти начальное приб.пижение к у с помощью явной формулы Эйлера

У° = Уп-1 +М(ж 1, y -i),

однако это не всегда целесообразно. В случае f(x, у) = Л/у имеем у° = (1 --M/i)y i и при \Mh\ 3> 1 может случиться, что такое приб.пижение будет слишком далеко от истинного решения. Поэтому более безопасный, но не всегда самый эффективный вариант -это положить у° = y i. Перепишем (13) в виде системы не.пииейных уравнений относительно у :

У - Уп-1 - /if (Жгг, Угг) = О,

и применим итерационный метод Ньютона. В данном конкретном случае интерполяционная формула Ньютона приобретает вид

У = Уп - - Щ{п, yi)) {yi - Уп-1 - hfixn, yfj), (14)

где к - номер итерации.

В одном из методов решения жестких систем за у принимается значение у получаемое из (14) при у° = Уп-\- Тогда имеем

Уп = yJ, = Уп-1 + hE - MyiXn, Уп)У{{хп, Уп-l).

Рассмотрим случай скалярного уравнения у = My; тогда Mh 1

Уп = Уп-1 + Z-Т7тУп-1

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 [132] 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208