Ващенко Г.В., Добронец Б.С. Надежность информационных систем. Лабораторный практикум

Подождите немного. Документ загружается.

Предположим, что f

(x) на интервале [a, b] не меняет знак и для опре-

деленности f

(x) > 0. Тогда k можно положить

k = max

[a,b]

(x).

Рассмотрим систему линейных алгебраических уравнений вида

x = Bx + b, (35)

где B — заданная матрица n-го порядка, b ∈ R

— заданный вектор.

Приведем условия, при которых отображение Bx + b является сжимаю-

щим.

ρ(Bx − By) = ρ(Bx − y) ≤ αρ(x − y),

другими словами, необходимо выполнение условия:

sup

x,y∈R

ρ(Bx − y)

ρ(x − y)

= sup

v6=0

ρ(Bv − 0)

ρ(v − 0)

≤ α,

или

sup

v6=0

||(Bv||

||v||

≤ α . (36)

В пространстве квадратных матриц мы можем задать норму следующим об-

разом:

||B|| = sup

v6=0

||(Bv||

||v||

. (37)

Можно показать, что

||A|| ≤ (

i,j=1

)

1/2

. (38)

Следовательно, для сходимости метода простой итерации достаточно пока-

зать, что норма матрицы ||B|| удовлетворяет следующему неравенству ||B|| ≤

α или выполнено более сильное условие

(

i,j=1

)

1/2

≤ α.

Системы нелинейных уравнений

Рассмотрим систему нелинейных уравнений с

неизвестными:

= ϕ

, x

, . . . , x

= ϕ

, x

, . . . , x

), (39)

. . . . . .

= ϕ

, x

, . . . , x

или в векторном виде

x = ϕ(x).

Если отображение ϕ сжимающее, то для нахождения решения системы

(39) мы можем применить метод простых итераций

= ϕ(x

), x

= ϕ(x

), . . . , x

= ϕ(x

n−1

). (40)

Предположим, что вектор-функция ϕ имеет непрерывные частные производ-

ные по x

, x

, . . . , x

, и обозначим A = (a

), где

= max

Ω

|∂ϕ

/∂x

|. (41)

Пусть x, y ∈ Ω. Согласно формуле конечных приращений Лагранжа,

имеем

(x) − ϕ

(y) =

j=1

∂ϕ

(ξ

)

∂x

− y

Следовательно,

ϕ(x) − ϕ(y) = A(x − y)

ρ(ϕ

(x) − ϕ

(y)) = ||ϕ

(x − ϕ

(y)|| = ||A(x − y)|| ≤

≤ ||A||||x − y|| = ||A||ρ(x − y),

т. е.

ρ(ϕ

(x) − ϕ

(y)) ≤ ||A||ρ(x − y).

Если ||A|| ≤ 1, то оператор ϕ является сжимающим и итерационный процесс

(40) сходится.

В интервальном виде метод простой итерации записывается как

i=1

= ϕ(x

) ∩ x

, i = 0, . . .

Причем для начального приближения x

должно выполняться включение

∗

⊆ x

Метод Ньютона

Пусть f ∈ C

[a, b] — некоторая дважды непрерывно дифференцируемая

функция. Рассмотрим задачу нахождения корня уравнения:

f(x) = 0. (42)

Теорема Пусть на отрезке [a, b] производная f

(x) не меняет знак

и, кроме того, f(a)f(b) < 0, тогда на отрезке существует корень x

∗

уравнения (42), причем единственный.

Доказательство этой теоремы очевидно из геометрических соображений.

Приведем следующую цепочку тождественных преобразований:

−f(x

) = f(x

∗

) − f(x

) = f

(ξ)(x

∗

− x

), x

∈ [a, b].

Разделив обе части равенства на f

(ξ), получаем

∗

− x

= −f(x

)/f

(ξ),

или

∗

= x

− f(x

)/f

(ξ), ξ ∈ [x

∗

, x

Таким образом, если точка x

близка к x

∗

, то в последнем равенстве можно

заменить неизвестное значение ξ на x

. Получаем приближенное равенство

∗

≈ x

− f(x

)/f

На основе этого приближения мы можем построить следующий итерацион-

ный процесс:

= x

i−1

− f(x

i−1

)/f

i−1

), i = 1, 2, . . . , (43)

который называют методом Ньютона или — поскольку (43) геометрически

означает уравнение касательной к графику функции f в точке x

i−1

— мето-

дом касательных.

Оценим скорость сходимости. Согласно формуле Тейлора имеем

0 = f(x

∗

) = f(x

i−1

) − f

i−1

)(x

∗

− x

i−1

) + f

(ξ)

∗

− x

i−1

)

или

∗

= x

i−1

− f(x

i−1

)/f

i−1

) + (x

∗

− x

i−1

)

(ξ)

i−1

)

, ξ ∈ [x

∗

, x

i−1

Вычитая из последнего равенства (43) и предполагая

(ξ)

i−1

)

≤ β,

получаем

∗

− x

| ≤ β(x

∗

− x

i−1

)

Отсюда следует, что если выполнено неравенство

β|x

∗

− x

| < 1, (44)

то погрешность убывает очень быстро по квадратичному закону. Через n

итераций будем иметь

∗

− x

i+n

| ≤

(β(x

∗

− x

))

В качестве примера рассмотрим следующее уравнение:

f(x) = x

− a = 0, f

(x) = 2x.

Метод Ньютона запишется в виде

i+1

= x

−

− a

− a/x

), (45)

за начальное приближение можно взять точку x

= a. Тем самым мы по-

лучили известную формулу для вычисления квадратных корней. Заметим,

что условие (44) может не выполняться, однако мы можем интерпретировать

формулу (45) как метод простой итерации. Несложно убедиться, что условие

сходимости выполнено и итерационный процесс (45) будет сходиться. Внача-

ле скорость сходимости будет как у метода простой итерации. В дальнейшем

при выполнении условия (44) скорость сходимости будет квадратичная и про-

цесс очень быстро сойдется.

Перейдем к интервальным методам Ньютона. Пусть f ∈ C

и известно,

что x

∗

∈ x

— корень уравнения f(x

∗

) = 0. Пусть f

(ξ) 6= 0, ∀ξ ∈ x

Интервальный метод Ньютона можно записать в виде

mid x

k+1

= (mid x

− f(mid x

)/f

)) ∩ x

. (46)

Последовательность {x

}, вычисленная по формулам (46), обладает свой-

ствами:

⊃ x

, . . . ,

lim

k→∞

= x

∗

, x

∗

∈ x

, ∀k.

Метод Кравчика

Предположим, что функция f удовлетворяет условиям метода Ньютона.

Рассмотрим оператор K:

K(x, x, y) = x − yf(x) + (1 − yf

(x))(x − x).

Теорема Пусть x

∗

∈ x

. Определим последовательность {x

}

k+1

= K(x

, x

, y

) ∩ x

где x

∈ x

, y

произвольны. Тогда:

1) x

∗

∈ x

, ∀k;

2) wid x

→ 0 при условии

1 − yf

) ⊃ [0, q], q < 1, ∀k.

На практике в качестве x

часто берут mid x

, а y

≈ (mid f

))

−1

или

(mid f

))

−1

. Заметим, что метод Кравчика не требует обращения интер-

вальных матриц.

Интерполяционный многочлен Лагранжа

Пусть известны значения некоторой функции f, заданной на отрез-

ке [a, b] в n + 1 различных точках a = x

< x

... < x

= b, которые

обозначим следующим образом:

= f(x

), i = 0, 1, ..., n.

Множество точек {x

} далее будем называть сеткой. Возникает задача при-

ближенного восстановления функции f в произвольной точке x. Построим

для этого алгебраический многочлен L

(x) степени n:

) = f(x

), i = 0, 1, ..., n. (47)

Такой многочлен будем называть интерполяционным , точки x

узлами

интерполяции. Рассмотрим общий вид интерполяционного многочлена:

(x) =

i=0

(x)f

, (48)

где

(x − x

)...(x − x

i−1

)(x − x

i+1

)...(x − x

)

− x

)...(x

− x

i−1

)(x

− x

i+1

)...(x

− x

)

, i = 0, 1, ..., n. (49)

Приведем примеры интерполяционного многочлена Лагранжа. Пусть n =

1, тогда

(x) =

x − x

− x

x − x

− x

. (50)

При n = 2

(x) =

(x − x

)(x − x

)

− x

)(x

− x

)

(x − x

)(x − x

)

− x

)(x

− x

)

(x − x

)(x − x

)

− x

)(x

− x

)

Погрешность интерполяции.

Оценим функцию R

(x) = f(x) − L

(x). (51)

Предположим, что f ∈ C

n+1

, т.е. n + 1 раз непрерывно дифференцируема.

Представим R

(x) в следующем виде:

(x) = ω

(x)r

(x), (52)

где

(x) = (x − x

)(x − x

) . . . (x − x

). (53)

Зафиксируем произвольное значение x, x 6= x

и рассмотрим следующую

функцию от t :

ϕ(t) = L

(t) + ω

(t)r

(x) − f(t). (54)

Эта функция обращается в нуль при t = x

, t = x, т.е. в n + 2 точках.

По теореме Ролля производная по t от ϕ обращается в нуль по крайней

мере в n + 1 точке, вторая производная равна нулю в n точках и т.д. Таким

образом, существует точка ξ такая, что ϕ

(n+1)

(ξ) = 0 . Поскольку L

(n+1)

(ξ) =

0, ω

(n+1)

(ξ) = (n + 1)!, получаем

(n + 1)!r

(x) − f

(n+1)

(ξ) = 0.

Следовательно,

(x) =

(n+1)

(ξ)

(n + 1)!

Остаточный член выглядит следующим образом:

(x) = ω

(x)

(n+1)

(ξ)

(n + 1)!

Если известна оценка f

(n+1)

(ξ) ≤ M

n+1

, то

|f(x) − L

(x)| ≤ |ω

(x)|

(n+1)

(n + 1)!

. (55)

В качестве примера рассмотрим оценку погрешности при линейной ин-

терполяции n = 1. Заметим, что ω

(x) = (x − x

)(x − x

) и, следовательно,

max

]

|ω

(x)| = max

]

|(x − x

)(x − x

)| = h

/4, (56)

где h = x

− x

. Таким образом,

max

]

|f(x) − L

(x)| ≤ h

. (57)

0.1 Квадратурные формулы

Для нахождения численных значений интегралов вида

I =

f(x)dx (58)

введем понятие квадратурной формулы.

Определение. Приближенное равенство

f(x)dx ≈

j=1

f(x

), (59)

где q

– некоторые числа, называемые весами , x

- некоторые точки отрезка

[a, b], называется квадратурной формулой.

Говорят, что квадратурная формула точна для многочленов степени

m, если при замене f на произвольный алгебраический многочлен степени m

приближенное равенство (59) становится точным.

Одним из способов построения квадратурных формул является замена

подынтегральной функции f ее некоторым приближением или аппроксима-

цией, от которой достаточно легко вычислить интеграл. Чем точнее мы най-

дем аппроксимацию, тем точнее мы вычислим интеграл.

Рассмотрим наиболее простые квадратурные формулы.

Формула прямоугольников

В этой формуле мы аппроксимируем нашу функцию f постоянным зна-

чением, взятым в середине интервала:

f(x)dx ≈ (b − a)f((b − a)/2). (60)

Интегрируя разложение функции f в ряд Тейлора, несложно убедиться, что

f(x)dx = (b − a)f((b − a)/2) +

(b − a)

(ξ), ξ ∈ [a, b], (61)

т.е. мы получили формулу прямоугольников с остаточным членом .

Формула трапеций

Для построения этой формулы мы заменяем функцию f полиномом

Лагранжа первой степени, или, другими словами, линейной интерполяцией

f(x)dx ≈ (b − a)

f(a) + f(b)

. (62)

Остаточный член для формулы трапеций записывается в следующем виде:

(b − a)

/12f

(ξ).

Формула Симпсона

Следуя рассмотренной выше схеме построения квадратурных формул,

перейдем от линейной интерполяции к интерполяционному многочлену Лагран-

жа второй степени. Как было показано выше, для его построения необходимо

знать значения функции в трех точках. Пусть это точки x

−1

, x

и значе-

ния функции соответственно f

−1

, f

, далее для простоты предположим

h = x

− x

−1

= x

− x

. Тогда несложно убедиться, что интерполяционный

многочлен Лагранжа второй степени задается уравнением

y = f

− f

x +

−1

− 2f

+ f

−1

Отсюда легко находим

−h

ydx = h(f

−1

+ 4f

+ f

)/3.

Таким образом, формула Симпсона, называемая также формулой парабол,

имеет вид

−h

f(x)dx ≈ h(f

−1

+ 4f

+ f

)/3. (63)

Остаточный член для формулы Симпсона записывается в следующем виде:

/90f

(4)

(ξ).

Усложненные квадратурные формулы

Как видно из приведенных остаточных членов, при больших интервалах

[a, b] точность квадратурных формул может быть низкой. Для повышения

точности применяется следующий прием: интервал интегрирования разбива-

ется на подинтервалы и на каждом подинтервале применяется своя квадра-

турная формула. Остановимся подробнее на применении формулы трапеций.

Разобьем интервал интегрирования на подинтервалы [x

, x

i+1

] таким образом,

что

a = x

< x

... < x

= b, h = x

i+1

− x

Имеем

i+1

f(x)dx = h(f

+ f

i+1

)/2 + h

/12f

(ξ). (64)

Суммируя по всем подинтервалам, получаем усложненную формулу трапе-

ций

f(x)dx ≈ h(f

/2 + f

+ f

+ . . . + f

n−1

+ f

/2) = I

(65)

с остаточным членом h

(b−a)

(ξ). Аналогично получается и усложненная

формула Симпсона

f(x)dx ≈ h(f

+ 4f

+ 2f

+ 4f

. . . + 4f

2n−1

+ f

)/3 = I

(66)

с остаточным членом h

(b−a)

180

(ξ).

Из выражений остаточных членов для квадратурных формул видно, что

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

пени, формула Симпсона — для многочленов третьей степени.

Погрешность формулы трапеций и формулы Симпсона удовлетворяет

неравенствам

|I − I

| ≤ h

(b − a)

max

[a,b]

(x)|,

|I − I

| ≤ h

(b − a)

180

max

[a,b]

(x)|.

Методы решения обыкновенных дифференциальных уравне-

ний

Задача Коши для обыкновенных дифференциальных уравне-

ний

Рассмотрим задачу Коши для обыкновенного дифференциального

уравнения первого порядка

= f(x, y), x ∈ (0, l), y(0) = y

, (67)

т.е. найдем решение уравнения (67), удовлетворяющее начальному условию

y(0) = y

. В дальнейшем будем считать, что решение задачи (67) существу-

ет и единственно. Нахождение точного аналитического решения задачи (67)

в общем случае не представляется возможным, поэтому мы рассмотрим на-

хождение приближенного численного решения на y

на сетке

= {0 = x

< x

... < x

= l}, h = x

− x

i−1

Далее приближенное решение y

) будем обозначать y

, положим y

= y

Численное решение будем искать последовательно. Зная значение функ-

ции y

в точке x

, найдем приближенное численное значение y

i+1

в точке

i+1

. Обратим внимание на тот факт, что если в точке x

мы знаем значе-

ние y

, то, следовательно, можем вычислить производную y

) = f(x

, y

Воспользуемся формулой Тейлора

y(x

) = y(x

) + y

)(x

− x

) + y

(ξ)(x

− x

)

/2, ξ ∈ (x

− x

или

y(x

) = y(x

) + hf(x

, y

) + y

(ξ)h

/2. (68)

Если h достаточно мало, то неизвестным слагаемым y

(ξ)h

/2 можно прене-

бречь:

y(x

) ≈ y(x

) + hf(x

, y(x

)).

Метод Эйлера

Приближенное решение задачи Коши будем вычислять на сетке ω

по

формуле

= y

i−1

+ hf(x

i−1

, y

i−1

). (69)

Заметим, что на каждом шаге мы ошибаемся на величину порядка y

(ξ)h

/2.

Следует ожидать, что глобальная ошибка |y(x

) − y

| будет ограничена сум-

мой локальных ошибок y

(ξ)h

/2. Более тонкая оценка приводится в лемме.

Лемма 2. Пусть M

= max |f(x, y)|, M

= max |f

(x, y)|, M

= max |f

(x, y)|

, тогда

|y(x

) − y

| ≤ h

− 1).

В качестве примера рассмотрим простейшее дифференциальное уравне-

ние

= y, y

= 1,

т.е. f(x, y) = y и точное решение этого уравнения y(x) = e

. Тогда

= y

+ hy

= 1 + h,

= y

+ hy

= (1 + h)y

= (1 + h)

Аналогично не сложно убедиться, что y

= y

k−1

+ hy

k−1

= (1 + h)y

k−1

(1 + h)

. Таким образом, полагая h = 1/n,

(1) = (1 + h)

= (1 +

)

Итак, мы получили известный предел для определения числа e ≈ (1 +

)

Одномерное параболическое уравнение

Обратимся к одномерному параболическому уравнению теплопроводно-

сти с коэффициентом, имеющим разрыв первого рода. Его можно решить

методом, основанным на вычислении невязки от специального сплайна, ап-

проксимирующего разностное решение. Сплайн построен так, чтобы на линии

разрыва коэффициента уравнения при фиксированном x выполнялись усло-

вия сопряжения (непрерывность решения и потока).

Рассмотрим уравнение

∂

u = ∂

(p∂

u) − qu + f, (70)

p(x, t) ≥ c

> 0, q(x, t) ≥ 0,

где (x, t) ∈ Q = (0, 1) ×(0, T ). Определим для функции u начальные и крае-

вые условия:

u(x, 0) = 0, x ∈ [0, 1], (71)