Ващенко Г.В., Добронец Б.С. Надежность информационных систем. Лабораторный практикум
Подождите немного. Документ загружается.
Пример 3. Пусть дана функция f = 2xy + z
2
+ xz −2x + 2z. Требуется
найти ее интервальное расширение на интервале x = [−1, 1] × [0, 1] × [0, 1].
Как и в примере 2, находим интервальные расширения частных производ-
ных f
0
x
= [−2, 1], f
0
y
= [−2, 2], f
0
z
= [1, 5]. Тогда x
∗
= [−1, 1] × [0, 1] × 1.
Повторно вычисляем f
0
x
= [−1, 1], f
0
y
= [−2, 2]. Сужение интервалов дальше
не происходит.
Переменная x удовлетворяет условиям алгоритма: не содержит квадрата
и wid(x
∗
) = 2. Приводим подобные члены по x. Получаем функцию f =
x(2y + z −2) + z
2
+ 2z, естественное интервальное расширение этой функции
на интервале x
∗
= [−1, 1] × [0, 1] × 1 совпадает с объединенным, поскольку
все интервальные переменные встречаются только один раз и только в первой
степени (z не интервальная переменная).
Заметим, что алгоритм построения интервального расширения можно
распространить и на более общие случаи полиномов степени m от n перемен-
ных. Покажем как это можно сделать на примере полинома третьей степени
от трех переменных
f(x) =
3
X
i,j,k=0
a
ijk
x
i
x
j
x
k
.
В п.1 алгоритма объединенное интервальное расширение уже не будет сов-
падать с их естественным интервальным расширением. Вследствие этого на-
личие монотонности по некоторым переменным может оказаться невыявлен-
ным. Для уточнения значений производных заметим, что они — полиномы
второй степени и необходимо установить либо отрицательность значения мак-
симума частной производной на интервале x, либо положительность значе-
ния ее минимума на интервале x. Ставя задачу таким образом, мы сводим ее
к уже упомянутому случаю — нахождению интервального расширения для
полиномов второй степени.
Таким образом, можно рекурсивно понижать степень полиномов и доби-
ваться значительного сужения ширины интервальных расширений.
Описанный выше подход можно распространить на случай произволь-
ных вещественных функций, для которых известны интервальные расшире-
ния g
i
. Предположим, что известны индексы i, при которых
0 6∈ g
i
(x). (10)
Это означает, что по этим переменным функция монотонна и, следовательно,
нижняя и верхняя границы функции достигаются при граничных значениях
интервала. Исходя из этого мы можем положить f(x) = f(z
1
), f(x) = f(z
2
),
51
где
z
1,i
=
x
i
, если g
i
> 0;
x
i
, если g
i
< 0;
x
i
, если g
i
3 0,
z
2,i
=
x
i
, если g
i
< 0;
x
i
, если g
i
> 0;
x
i
, если g
i
3 0.
Таким образом, вычисление каждой из границ интервального расширения
сводится к вычислению интервальных расширений той же функции, но при
более узких значениях интервального аргумента. Далее к каждому такому
вычислению интервального расширения мы можем вновь применить одну из
перечисленных выше процедур.
Пример 4. Рассмотрим следующую функцию трех переменных, задан-
ную на [0, 1]
3
, a = b = c = 0.5:
f(x
1
, x
2
, x
3
) = ax
1
x
2
x
3
+ 2x
1
x
2
+ 2bx
2
x
3
+ cx
1
− (1 − b)x
2
− bx
3
.
Несложно подсчитать, что
∂
1
f(x
1
, x
2
, x
3
) = [c, 2 + a + c] > 0,
∂
2
f(x
1
, x
2
, x
3
) = [−1 − b, 1 + a + b] 3 0,
∂
3
f(x
1
, x
2
, x
3
) = [−b, a + b] 3 0.
Поскольку f(x
1
, x
2
, x
3
) монотонно возрастает по x
1
, то для нахождения f мы
можем положить x
1
= x
1
= 1. Относительно поведения функции по второму
и третьему аргументам мы ничего определенного сказать не можем. В этом
случае
∂
2
f(x
1
, x
2
, x
3
) = [1 − b, 1 + a + b] > 0,
∂
3
f(x
1
, x
2
, x
3
) = [−b, a + b] 3 0.
Далее при x
1
= x
1
= 1 в силу положительности частной производной по
второму аргументу f(x
1
, x
2
, x
3
) возрастает по x
2
, следовательно, можно по-
ложить x
2
= x
2
= 1. При фиксированных значениях x
1
= 1, x
2
= 1
∂
3
f(x
1
, x
2
, x
3
) = a + b > 0.
Таким образом, нам удалось последовательно построить верхнюю границу
интервального расширения в виде
f(x
1
, x
2
, x
3
) = f(x
1
, x
2
, x
3
).
Аналогично можно построить и нижнюю границу
f(x
1
, x
2
, x
3
) = f(x
1
, x
2
, x
3
).
Несложно видеть, что построенное интервальное расширение совпадает с объ-
единенным интервальным расширением.
52
Прямые методы
Рассмотрим систему линейных алгебраических уравнений
Ax = b, (11)
где b = (b
i
), i = 1, . . . , n — известный вектор, A = (a
i,j
), i, j = 1, . . . , n —
невырожденная матрица. Тогда вектор x = A
−1
b — решение системы (11).
Предположим, что A и b содержат ошибки и известно, что их элементы при-
надлежат соответствующим интервальным числам
a
i,j
∈ a
i,j
,
b
i
∈ b
i
, i, j = 1, . . . , n.
Пусть также, все матрицы из множества A не вырождены. Множество
векторов
X = {x|Ax = b, A ∈ A, b ∈ b}
назовем множеством решений системы
Ax = b. (12)
Пример 5. Пусть необходимо решить систему интервальных линейных
алгебраических уравнений
Ax = b
с матрицей A и правой частью b.
A =
[2, 4] [−2, 1]
[−1, 2] [2, 4]
, b =
[−2, 2]
[−2, 2]
. (13)
Минимальным интервальным вектором, содержащим множество ее ре-
шений, является интервальный вектор x = ([−4, 4], [−4, 4])
T
.
Множество X может быть описано следующим образом:
X = {x|x ∈ R
n
, Ax ∩ b 6= ∅} (14)
или
{x|x ∈ R
n
, 0 ∈ Ax − b}.
Справедливо следующее утверждение.
Замечание.
x ∈ X ⇔ |mid Ax − mid b| ≤ rad(b).
53
Поставим задачу найти интервальный вектор x ∈ R
n
, содержащий мно-
жество X.
Самым простым способом для нашего примера будет метод Крамера.
Действительно, для СЛАУ A = (a
ij
), b = (b
i
), i, j = 1, 2. Решение находится
по формулам
x
1
=
∆
1
∆
, x
2
=
∆
2
∆
,
где
∆ = a
11
a
22
− a
12
a
21
,
∆
1
= b
1
a
22
− b
2
a
21
,
∆
2
= a
11
b
2
− a
12
b
1
.
Для решения СЛАУ с интервальными коэффициентами из примера (13) сде-
лаем интервальное расширение формул Крамера. Непосредственными вычис-
лениями получаем x
1
= [−6, 6], x
2
= [−6, 6].
Как видим, ответ несколько шире оптимального. Этот факт легко объ-
ясняется, поскольку рациональные выражения для x
1
, x
2
не удовлетворяют
условию теоремы , т.е. содержат переменные более одного раза.
Применять формулы Крамера для решения больших систем крайне не
рационально. В вычислительной математике для этих целей используют пря-
мые методы: метод Гаусса, LU- и QR-разложения и т. д.
Рассмотрим еще ряд примеров.
Пример 6. Пусть необходимо решить систему интервальных линейных
алгебраических уравнений
A =
[2, 4] [−1, 1]
[−1, 1] [2, 4]
, b =
[−3, 3]
[0, 0]
. (15)
Множество векторов X этой задачи согласно замечанию 1.
X = {x|x ∈ R
2
, 2|x
2
| ≤ |x
1
|, 2|x
1
| ≤ 3 + |x
2
|}
Пример 7. Пусть
A =
2 α
β 2
, b =
1.2
−1.2
, α, β ∈ [0, 1]. (16)
Заметим, что по правилу Крамера
x
1
= 1.2(2 − α)/(4 − αβ),
x
2
= −1.2(2 − β)/(4 − α β).
54
Следовательно,
2X =
[0.3, 0.6]
[−0.6, −0.3]
.
Применяя формально правило Крамера к системе (20), мы получаем интер-
вальный вектор ([0.3, 1.2], [−0.8, 1.2]).
Найдем обратную матрицу
A
−1
= 2
1
4 − αβ
2 α
β 2
¯
¯
¯
¯
¯
α, β ∈ [0, 1]
=
=
[1/2, 2/3] [0, 1/3]
[0, 1/3] [1/2, 2/3]
.
Несложно убедиться, что
A
−1
b 6= 2X.
Пусть A — симметричная матрица
A =
2 α
α 2
, α ∈ [0, 1].
Множество решений системы линейных алгебраических уравнений с симмет-
ричной матрицей представляет собой отрезок с концами (0.4, −0.4), (0.6, −0.6).
Интервальная оболочка множества решений с симметричной матрицей
2X
sym
=
[0.4, 0.6]
[−0.6, −0.4]
6= 2X.
Метод Гаусса и LU-разложение
Интервальный метод Гаусса представляет корректное интервальное рас-
ширение метода Гаусса для СЛАУ.
Рассмотрим систему интервальных линейных алгебраических уравнений
с квадратной матрицей размерности n:
a
11
, a
12
, . . . a
1n
a
21
, a
22
, . . . a
2n
.
.
.
.
.
.
.
.
.
a
n1
, a
n2
. . . a
nn
x
1
x
2
.
.
.
x
n
=
b
1
b
2
.
.
.
b
n
. (17)
55
Метод Гаусса состоит из двух ходов: прямого и обратного. На прямом ходе с
помощью элементарных преобразований будем последовательно стремиться
привести систему (17) к верхнетреугольному виду.
Первый шаг. Преобразуем систему к виду
a
11
, a
12
, . . . a
1n
0, a
(1)
22
, . . . a
(1)
2n
.
.
.
.
.
.
.
.
.
0, a
(1)
n2
. . . a
(1)
nn
x
1
x
2
.
.
.
x
n
=
b
1
b
(1)
2
.
.
.
b
(1)
n
, (18)
где
a
(1)
ik
= a
ik
− l
i1
a
1k
,
b
(1)
i
= b
i
− l
i1
b
1
, l
i1
= a
i1
/a
11
.
Далее переходим к подсистеме с матрицей A
(1)
и правой частью b
(1)
раз-
мерности n − 1. Применим к этой системе первый шаг и т. д. Окончательно
приводим cиcтему (17) к верхнетреугольному виду
a
11
a
12
. . . a
1n
0 a
(∗)
22
. . . a
(∗)
2n
.
.
.
.
.
.
.
.
.
0, 0 . . . a
(∗)
nn
x
1
x
2
.
.
.
x
n
=
b
1
b
(∗)
2
.
.
.
b
(∗)
n
. (19)
Обратный ход. Далее последовательно вычисляем неизвестные x
n
, x
n−1
,...,x
1
:
x
n
= b
(∗)
n
/a
(∗)
nn
,
x
i
=
b
i
−
P
i+1
j=n
a
(∗)
j
x
j
a
(∗)
ii
i = n − 1, ..., 1.
Интервальный метод Гаусса не всегда можно реализовать даже для регуляр-
ных матриц A.
Пример Райхмана
Рассмотрим интервальную матрицу
S(a) =
1 [0, a] [0, a]
[0, a] 1 [0, a]
[0, a] [0, a] 1
.
Проводя прямой ход метода Гаусса, получаем
S(a) ∼
1 [0, a] [0, a]
0 1 [−a
2
/(1 − a
2
), a/(1 − a
2
)]
0 0 [1 − a
2
− a
2
/(1 − a
2
), 1 + a
3
(1 − a
2
)]
.
56
Несложно убедиться, для любого a ∈ [(
√
5 − 1)/2, 1] следует, что 0 ∈ [1 −
a
2
− a
2
/(1 − a
2
), 1 + a
3
(1 − a
2
)]. Таким образом, метод Гаусса не может быть
закончен.
Терема 3. Если матрица A — M-матрица или имеет диагональное
преобладание, то решение x
G
по методу Гаусса существует и
2X ⊆ x
G
.
Согласно методу Холецкого матрицу A можно представить в виде про-
изведения двух треугольных матриц L, U .
L =
l
11
, 0 . . . 0
l
21
, l
22
. . . 0
.
.
.
.
.
.
l
n1
, l
n2
. . . l
nN
, U =
1, u
12
. . . u
1n
0, 1 . . . u
2n
.
.
.
.
.
.
0, 0 . . . 1
.
Коэффициенты матриц L, U находят по формулам:
l
i1
= a
i1
для i = 1, ..., n,
u
1i
= a
1i
/l
11
для i = 2, ..., n.
Далее справедливы рекуррентные соотношения:
l
is
= (a
is
−
s−1
X
k=1
l
ik
u
ki
), i = s, ..., n; s = 2, 3, ..., n,
u
si
= (a
is
−
s−1
X
k=1
l
ik
u
ki
)/l
ss
, i = s + 1, ..., n; s = 2, ..., n − 1.
Если матрицы L, U построены, предлагается следующий алгоритм. Сна-
чала решают вспомогательную систему
Ly = b.
Так как матрица L треугольная, то нетрудно выписать рекуррентные
соотношения
y
1
= b
1
/l
11
, y
2
= (b
2
− l
12
y
1
)/l
22
, ...,
y
n
= (b
n
−
n−1
X
i=1
l
ni
y
i
)/l
nn
.
Зная вектор y, находят решение системы (17) с помощью треугольной мат-
рицы U.
Ux = y.
57
Аналогично выписывают рекуррентные соотношения
x
n
= y
n
, x
n−1
= y
n−1
− u
n−1,n
x
n
, ...,
x
1
= y
1
−
n
X
i=2
u
1i
x
i
.
Метод LU-разложения наиболее предпочтителен, если систему (17) необхо-
димо решать для разных правых частей b. Один раз нужно найти матрицы
L, U, а затем пользоваться обратным ходом нахождения векторов Y, X.
Пример 8. Рассмотрим
A =
2 [−1, 0]
[−1, 0] 2
, b =
1.2
−1.2
. (20)
A — M-матрица с LU-разложением
L =
1 0
[−0.5, 0] 1
, U =
2 [−1, 0]
0 [1.5, 2]
.
Следовательно,
x =
[0.2, 0.6]
[−0.8, −0.3]
.
Метод простой итерации
Пусть в некоторой области Ω ⊂ R
n
задано отображение
ϕ : R
n
→ R
n
, причем для любых x, y ∈ Ω выполнено неравенство
ρ(ϕ(x) − ϕ(y)) ≤ αρ(x − y), (21)
где 0 ≤ α < 1. Отображение, удовлетворяющее свойству (21), называет-
ся сжимающим отображением. Точка x
∗
называется неподвижной точкой
отображения ϕ, ϕ(x
∗
) = x
∗
. Другими словами, неподвижная точка — это
решение уравнения
ϕ(x) = x.
Теорема.[о принципе сжимающих отображений] Всякое сжимающее отоб-
ражение имеет одну и только одну неподвижную точку x
∗
, причем
x
∗
= lim
k→∞
x
k
, (22)
где x
0
∈ Ω — произвольное начальное приближение,
x
k
= ϕ(x
k−1
). (23)
58
Доказательство. Покажем, что последовательность {x
k
} фундаменталь-
ная. Действительно, рассмотрим процесс, который называется методом про-
стой итерации:
x
1
= ϕ(x
0
), x
2
= ϕ(x
1
) = ϕ
2
(x
0
); x
n
= ϕ(x
n−1
) = ϕ
n
(x
0
).
Будем считать для определенности n ≤ m, имеем
ρ(x
n
− x
m
) = ρ(ϕ
n
(x) − ϕ
m
(x)) ≤ α
n
ρ(x
0
− x
m−n
) ≤
≤ α
n
(ρ(x
0
− x
1
) + ρ(x
1
− x
2
) + . . . + ρ(x
m−n−1
− x
m−n
)) ≤
≤ α
n
ρ(x
0
− x
1
)(1 + α + α
2
+ . . . + α
m−n−1
) ≤ α
n
ρ(x
0
− x
1
)/(1 − α).
Поскольку α < 1, то при достаточно большом n величина ρ(x
n
− x
m
) сколь
угодно мала, то последовательность {x
k
} фундаментальная и предел (22)
существует. Тогда в силу непрерывности ϕ
ϕ(x
∗
) = ϕ( lim
k→∞
x
k
) = lim
k→∞
ϕ(x
k
) = lim
k→∞
x
k+1
= x
∗
.
Итак, существование неподвижной точки доказано, покажем ее единствен-
ность. Если
ϕ(x) = x, ϕ(y) = y,
то неравенство (21) принимает вид
ρ(x − y) ≤ αρ(x − y).
Следовательно, ρ(x − y) = 0 , или x = y. Теорема доказана.
Рассмотрим простейшие применения принципа сжимающих отображе-
ний к решению нелинейных уравнений.
Пусть задано уравнение с одной неизвестной x:
x = ϕ(x), (24)
где ϕ(x) — заданная функция. Уравнение (24) может иметь различное число
корней или совсем не иметь решений.
Функция ϕ удовлетворяет на отрезке [a, b] условию Липшица с постоян-
ной α, если для любых x
1
, x
2
∈ [a, b] выполняется неравенство
|ϕ(x
1
) − ϕ(x
2
)| ≤ α|x
1
− x
2
|. (25)
Если функция ϕ(x) дифференцируема на отрезке [a, b], то она удовле-
творяет на [a, b] условию Липшица с постоянной
α ≤ max
[a,b]
|ϕ
0
(x)|. (26)
59
Теорема. Пусть функция ϕ удовлетворяет на отрезке [a, b] условию
Липшица с постоянной α, причем
0 ≤ α < 1. (27)
Тогда уравнение (24) имеет на отрезке [a, b] единственное решение
x
∗
= lim
k→∞
x
k
, (28)
где x
0
— начальное приближение из отрезка [a, b],
x
k
= ϕ(x
k−1
). (29)
Доказательство непосредственно вытекает из принципа сжимающих отоб-
ражений.
Оценим скорость сходимости (29). Имеем
x
∗
= ϕ(x
∗
). (30)
Вычтем (29) из (30), тогда
|x
∗
− x
k
| ≤ α
k
|x
∗
− x
0
| ≤ α
k
(b − a). (31)
Таким образом, мы видим, что расстояние между точным решением нашего
уравнения и приближенным убывает в геометрической прогрессии.
На практике уравнения вида (24) встречаются довольно редко. Обычно
нелинейные уравнения задаются в виде
f(x) = 0. (32)
Приведем способ преобразования (32) к виду (24). Заметим, что для лю-
бого k 6= 0 следующее уравнение эквивалентно предыдущему:
f(x)
k
= 0 (33)
и
x = x −
f(x)
k
. (34)
Следовательно, мы можем положить ϕ(x) = x −
f(x)
k
. В силу замечания 25
для построения сходящегося процесса нам необходимо выполнение условия
(27):
|ϕ
0
(x)| = |1 −
f
0
(x)
k
| ≤ 1.
60