Лабораторная работа по численному анализу: машинная точность, ряды, ОДУ

Санкт-Петербургский политехнический университет Петра Великого
Физико-механический институт
Высшая школа прикладной математики и вычислительной физики
Направление подготовки
"01.03.02. Прикладная математика и информатика"
Контроль надежности инженерных и научных вычислений
Лабораторная работа №1
Работу выполнил:
Хаяйнен Н.К.
Группа: 5030102/20101
Преподаватель:
к.ф.-м.н. Репин Сергей
Игоревич
Санкт-Петербург
2026
Содержание
1 Задача 1(a): Определить машинную точность своего ПК
1.1 Постановка задачи . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
1.2 Результат . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3
3
3
2 Задача 1(b): Чувствительность суммирования к критерию останова
2.1 Постановка задачи . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.2 Результаты . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.3 Выводы . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3
3
4
5
3 Задача 1(c): Составная задача (СЛАУ + ОДУ)
3.1 Постановка задачи . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.2 Аналитическое решение . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.2.1 Матрица системы . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.2.2 Вектор правой части . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.2.3 Решение ОДУ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.3 Численное решение . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.4 Выводы . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
6
6
6
6
6
7
7
7
4 Общие выводы
8
2
1
Задача 1(a): Определить машинную точность своего ПК
1.1
Постановка задачи
Из лекций (стр. 15):
«The largest real number 𝜀 > 0 such that 1 + 𝜀 = 1 is called the machine accuracy of the computer.»
Задача 1(a): Define the machine accuracy of your PC.
Перевод: определить машинную точность своего ПК — наибольшее вещественное 𝜀 > 0, такое
что 1 + 𝜀 = 1 в арифметике с плавающей точкой.
1.2
Результат
Машинная точность была найдена следующим циклом:
eps = 1.0
while (1.0 + eps) != 1.0:
eps /= 2.0
Полученное значение:
𝜀 = 1,110223 × 10−16
Это совпадает с теоретическим значением для 64-битной арифметики IEEE 754 (мантисса
52 бита, 𝜀теор = 2−52 ≈ 2,22 · 10−16 ).
2
Задача 1(b): Чувствительность суммирования к критерию останова
2.1
Постановка задачи
Из лекций (стр. 24):
Задача 1(b): Verify sensitivity of summation with respect to stopping criterion for a series having
low convergence rate.
Перевод: верифицировать чувствительность суммирования к критерию останова для ряда с
медленной скоростью сходимости.
Конкретно: вычислить знакопеременный ряд (функцию Дирихле)
𝜂(𝑠) ≈ 𝜂𝑁 (𝑠) :=
𝑁
∑︁
(−1)𝑘+1
𝑘=1
𝑘𝑠
используя критерий останова
|𝜂𝑁 (𝑠) − 𝜂𝑁 −1 (𝑠)| ≤ 𝜀,
что эквивалентно условию
1
≤ 𝜀.
𝑁𝑠
В работе рассмотрены следующие значения параметра:
1. 𝑠 = 1, точное значение 𝜂(1) = ln 2 ≈ 0,693;
3
,
2. 𝑠 = 1,1, эталонное значение 𝜂(1,1) ≈ 0,709;
3. 𝑠 = 1,5, эталонное значение 𝜂(1,5) ≈ 0,765;
4. 𝑠 = 2, эталонное значение 𝜂(2) ≈ 0,822.
2.2
точное понимаю, а эталонное нет
Грамотно писать так: точное значение
с точностью до 3 значащих цифр равно
Результаты
Таблица 1: Значения 𝜂𝑁 (𝑠) для альтернативного ряда при различных 𝑠 и 𝜀
𝑠
𝜀
𝑁
𝜂𝑁
|𝜂𝑁 − 𝜂ref |
ref ???
−1
1.0 10
10
0.6456349206349207 4.7512259925 · 10−2
1.0 10−2
100
0.6881721793101955 4.9750012497 · 10−3
1.0 10−3
1 000
0.6926474305598204 4.9975000012 · 10−4
−4
1.0 10
10 000
0.6930971830599454 4.9997500000 · 10−5
1.0 10−5
100 000
0.6931421805849454 4.9999749999 · 10−6
1.0 10−6 1 000 000 0.6931466805601952 4.9999975005 · 10−7
1.0 10−7 10 000 000 0.6931471305599485 4.9999996810 · 10−8
1.1 10−1
9
0.7506980515909494 4.1889211637 · 10−2
1.1 10−2
66
0.7038675883534551 4.9412516008 · 10−3
−3
1.1 10
534
0.7083097045170532 4.9913543722 · 10−4
1.1 10−4
4 329
0.7088588405455475 5.0000591270 · 10−5
−5
1.1 10
35 112
0.7088038500218552 4.9899324223 · 10−6
1.1 10−6
284 804
0.7088083499323503 4.9002192715 · 10−7
1.1 10−7 2 310 130 0.7088087999306065 4.0023670933 · 10−8
1.5 10−1
5
0.8033394182365932 3.8192393627 · 10−2
1.5 10−2
22
0.7604664927038867 4.6805319057 · 10−3
1.5 10−3
100
0.7646507743520380 4.9625025756 · 10−4
−4
1.5 10
465
0.7651968085966251 4.9783987030 · 10−5
1.5 10−5
2 155
0.7651520209186463 4.9963090512 · 10−6
−6
1.5 10
10 000
0.7651465246629077 4.9994668749 · 10−7
1.5 10−7
46 416
0.7651469746263961 4.9983199024 · 10−8
2.0 10−1
4
0.7986111111111112 2.3855922313 · 10−2
−2
2.0 10
10
0.8179621756109852 4.5048578131 · 10−3
2.0 10−3
32
0.8219939961054309 4.7303731867 · 10−4
2.0 10−4
100
0.8224175333741284 4.9500049970 · 10−5
−5
2.0 10
317
0.8224719933971315 4.9599730330 · 10−6
2.0 10−6
1 000
0.8224665339241130 4.9949998560 · 10−7
−7
2.0 10
3 163
0.8224670833854784 4.9961379811 · 10−8
Во всех примерах точность 10-7 по критерию оснановки дает точность
по решению 10-8. В чем разница в разных тестах?
4
Рис. 1: Зависимость ошибки |𝜂𝑁 − 𝜂ref | от критерия останова 𝜀 для различных значений 𝑘.
2.3
Выводы
При медленной сходимости ряда (при 𝑘, близком к 1) ошибка убывает достаточно медленно
при уменьшении 𝜀. Несмотря на то, что знакопеременный ряд сходится быстрее, чем соответствующий положительный ряд, для достижения высокой точности всё равно требуется
большое число слагаемых.
Например, при 𝑘 = 1 даже для 𝜀 = 10−7 требуется 𝑁 = 107 слагаемых, при этом ошибка
остаётся порядка 10−8 . Для значений 𝑘, близких к 1 (например, 𝑘 = 1,1), сходимость также
остаётся сравнительно медленной.
Это подтверждает, что критерий останова
|𝜂𝑁 − 𝜂𝑁 −1 | ≤ 𝜀
(то есть малость очередного члена ряда) не всегда гарантирует столь же быструю убывающую
погрешность, особенно при параметрах, близких к критическим значениям.
Я не вижу в этих примерах главного: одно и то же \eps в критерии
типа Рунге может приводить к совершенно разной точности по результату.
Сколько членов ряда надо тут не так важно.
5
3
Задача 1(c): Составная задача (СЛАУ + ОДУ)
3.1
Постановка задачи
Из лекций (стр. 26):
Задача 1(c): Функция 𝑧(𝑥) является решением дифференциального уравнения
𝑧 ′′ − 9𝑧 ′ − 10𝑧 = 0,
𝑧(0) = 1,
𝑧 = 𝑧(𝑥),
𝑥 ∈ [0, 3],
𝑁
𝑧 ′ (0) = − max |𝑎𝑖 |,
𝑖=1
где 𝑎 = {𝑎1 , . . . , 𝑎𝑁 } является решением системы 𝐵𝑎 = 𝑓 .
Коэффициенты матрицы 𝐵 ∈ R𝑁 ×𝑁 вычисляются как интегралы:
∫︁
)︁
2 𝜋 (︁
𝑏𝑖𝑗 =
sin(𝑖𝜉) sin(𝑗𝜉) + (𝑚 cos(𝑚𝜉) + 2 sin(𝑚𝜉))𝑒2𝜉 𝑑𝜉,
𝜋 0
𝑚 = 𝑖 + 𝑗 2.
Компоненты вектора 𝑓 :
𝑓𝑖 =
√︀
𝐻(𝑠𝑖 ),
𝜁(𝑠) =
𝐻(𝑠𝑖 ) =
∞
∑︁
1
𝑘=1
𝑘𝑠
,
ˆ 𝑖)
𝜁(𝑠
2𝑠𝑖
·
,
2𝑠𝑖 − 1 𝜁(𝑠𝑖 )
ˆ =
𝜁(𝑠)
𝑠𝑖 =
𝑖+1
,
𝑖
∞
∑︁
1
.
(2𝑘 + 1)𝑠
𝑘=0
Надо вычислить значение 𝑧(3) (для 𝑁 = 10, 20, 40), проинтегрировав уравнение численно
на некоторой сетке ℎ. Коэффициенты матрицы надо найти, применяя формулы численного
интегрирования и суммируя ряды численно.
3.2
Аналитическое решение
3.2.1
Матрица системы
В силу ортогональности системы {sin(𝑘𝜉)} на [0, 𝜋]:
∫︁
2 𝜋
sin(𝑖𝜉) sin(𝑗𝜉) 𝑑𝜉 = 𝛿𝑖𝑗 .
𝜋 0
Второй интеграл после двукратного интегрирования по частям:
∫︁
2
2 𝜋
(𝑚 cos(𝑚𝜉) + 2 sin(𝑚𝜉))𝑒2𝜉 𝑑𝜉 = 𝑒2𝜋 sin(𝜋𝑚) = 0,
𝜋 0
𝜋
Таким образом, 𝐵 = 𝐼 (единичная матрица).
3.2.2
Вектор правой части
ˆ
Используя соотношение 𝜁(𝑠)(1 − 2−𝑠 ) = 𝜁(𝑠):
𝐻(𝑠𝑖 ) =
2𝑠𝑖
2𝑠𝑖 − 1
·
= 1 =⇒ 𝑓𝑖 = 1.
2𝑠𝑖 − 1
2𝑠𝑖
Решение СЛАУ: 𝑎 = (1, . . . , 1)𝑇 , откуда 𝑧 ′ (0) = −1.
6
𝑚 ∈ N.
3.2.3
Решение ОДУ
Характеристическое уравнение: 𝛾 2 − 9𝛾 − 10 = 0, корни 𝛾1 = −1, 𝛾2 = 10.
Общее решение: 𝑧(𝑥) = 𝐶1 𝑒−𝑥 + 𝐶2 𝑒10𝑥 . Из НУ: 𝐶1 = 1, 𝐶2 = 0, т. е. 𝑧(𝑥) = 𝑒−𝑥 .
𝑧(3) = 𝑒−3 ≈ 0,04979.
3.3
Численное решение
Элементы 𝐵 вычислялись методом Симпсона (𝑛quad = 500), ОДУ решалось методом Эйлера и
Рунге-Кутта 4 (𝑛ODE = 1000).
В квадратурах есть регуляция точности: шаг
и в методе Рунге-Кутта, и в суммировании рядов есть
Таблица 2: Результаты
численного
решения
𝑧(3) на результат?
критерий
\eps. Как
они влияют
улучшает результат или нет?
𝑁 𝑧 ′ (0) числ. Их изменение/улучшение
Эйлер
РК-4 Точное
10
20
40
−1,024 −1,47 · 1010
−1,936 −5,85 · 1011
−346,6 −2,16 · 1014
−2,28 · 1010
−9,10 · 1011
−3,36 · 1014
0,04979
0,04979
0,04979
Рис. 2: Графики |𝑧(𝑥)| при различных 𝑁 . Красная линия — точное решение 𝑒−𝑥 , штриховые
— численные решения.
3.4
Выводы
Численное решение расходится с аналитическим на 10–14 порядков. Причина — плохая
обусловленность задачи: характеристические корни ОДУ 𝛾1 = −1 и 𝛾2 = 10 имеют разные
знаки, и малейшая погрешность 𝛿 в 𝑧 ′ (0) порождает ненулевую компоненту при 𝑒10𝑥 :
𝐶2 = −
𝛿
𝛿
=⇒ 𝑧(3) ≈ − 𝑒30 .
11
11
Даже при 𝛿 ∼ 0,02 ошибка достигает ∼ 1010 . Улучшение параметров численных методов не
может устранить эту принципиальную неустойчивость.
7
4
Общие выводы
1. Задача 1(a): машинная точность составляет 𝜀 ≈ 1,11 × 10−16 и определяет предел относительной точности вещественных вычислений.
2. Задача 1(b): для знакопеременного ряда с медленной сходимостью при 𝑘, близком к 1,
формальный критерий «следующий член мал» (1/𝑁 𝑘 ≤ 𝜀) не гарантирует столь же быстрого убывания погрешности суммирования; даже при 𝑘 = 1 для достижения точности
порядка 10−8 требуется ∼ 107 слагаемых.
3. Задача 1(c): на примере составной задачи показано, как малые погрешности при численном решении СЛАУ катастрофически усиливаются неустойчивым ОДУ (расхождение
на 10–14 порядков).
8