[ методы
и математическое
моделирование
Лекция 7
Численное решение:
алгоритмы, методы
и неприятности …
Ю.Н. Прошин ЧМММ. Лекция 7 #1
КУЛЬТУРА ВЫЧИСЛЕНИЙ НА ЭВМ
•До сих пор => Постановки задач и алгоритмы их решения.
•Однако, мы имеем цепочку «модель — алгоритм — программа».
•Одна из возможных причин несовпадения желаемого и
получаемого
=> несовпадение машинной арифметики с обычной из-за конечности разрядной сетки ЭВМ.
Возникающие ошибки могут привести к большим неприятностям, если их не контролировать и не соблюдать некоторые элементарные правила организации вычислений. Правила эти неформальны и напоминают правила хорошего тона. Уровень их выполнения определяет уровень вычислительной культуры пользователя ЭВМ. Поясним на примерах основные из этих правил.
Ю.Н. Прошин ЧМММ. Лекция 7 # Пример 1.
Коммутативность, ассоциативность, … Вычислим (на обычном калькуляторе или на ЭВМ) 1016 + 1 — 1016, 1032 + 1010 — 1016 — 1016 + 1, 1032 — 1032 + >> 10^16 + 1 - 10^16 >> single(10^8) + 1 - single(10^8) >> 10^16 - 10^16 + ans = ans = ans = 0 0 >> 10^32 + 10^10 - 10^32 >> 10^32 - 10^32 + 10^ ans = ans = 0 1.0000e+ >> >> Таким образом, в машинной арифметике нарушаются законы коммутативности и ассоциативности действий.
Применимость основных выводов элементарной математики ставится под сомнение.
Ю.Н. Прошин ЧМММ. Лекция 7 # Пример 2. Предел.
Известно, что >> n=single([1 1e5 1e7 1.1e7 1.2e7 2e7 3e7]) n= 1 100000 10000000 11000000 12000000 20000000 >> (1+1./n).^n ans = 2.0000 2.7220 3.2940 3.7110 4.1808 1.0000 1. Вывод: при вычислениях с ЭВМ применимость основного понятия высшей математики — предела — также ставится под сомнение.
Ю.Н. Прошин ЧМММ. Лекция 7 # Правило 1.
Определение. Машинным эпсилоном называется наименьшее представимое в ЭВМ число, удовлетворяющее условию Правило 1.
Величина М характеризует наименьшую относительную погрешность вычислений и зависит от конкретной ЭВМ и разрядности вычислений (single, double,... ).
Требовать БОЛЬШЕГО невозможно!
Ю.Н. Прошин ЧМММ. Лекция 7 # Задача и комментарий к правилу 1.
Задача. Найти минимальное число М, используемое компьютером при вычислениях по умолчанию, определить количество значащих цифр используемых при численных расчетах, выяснить возможности его увеличения (уменьшения). (На примере любого языка программирования Си, Паскаль, Фортран, Дельфи или вычислительного пакета MatLab, Maple, Mathematica, MathCad, Origin, Derive).
•Очевидно, если М > 10-k, то на данной ЭВМ нельзя гарантировать, что в результатах будет содержаться не менее k верных значащих цифр.
•Напомним, что цифра числа называется верной, если абсолютная погрешность числа не превосходит половины единицы того разряда, в котором эта цифра находится.
Чувствительность к исходным данным.
• Задача может быть чувствительна к малым ошибкам, допущенным при представлении исходных данных.
• Пример, корни уравнения p(x) = (x-2)2 = 0 равны двум, при изменении свободного члена на малую величину = 10- (x-2)2 = изменение в корнях много больше: x1,2 = 2 ± 10-3.
• Этот тип неустойчивости еще более выражен у полиномов более высокой степени. Корни следующего полиномиального уравнения p(x) = 0, где p(x) = (x-1)(x-2)...(x-20) = x20 - 210x19 + … суть реальные числа от 1 до 20 и хорошо разделены.
Задача: Измените коэффициент при x19: (-210) (-210 + 10-23) и проследите численно за катастрофическим изменением решения.
Пример 3. Эквивалентные формулы – Пусть М = 10-2 и требуется решить уравнение:
Округляем до двух значащих цифр.
Формулы для решения уравнения Верный ответ -- ошибка получена при вычитании близких Эквивалентная формула получаем Пример 4. Эквивалентные формулы – Оценка дисперсии случайной величины по измерениям :
Пусть х 1 = 12345.1, х 2 = 12345.2, х 3 = 12345.3.
При выборе формулы и порядка вычислений избегать вычитания близких чисел и деления на малые величины.
Столь большие отличия в ответах возникли из-за того, что матрицы коэффициентов систем 1 и 2 плохо обусловлены: определители i = det|ai| малы. Действительно, 1 = 0.001, 2 = -0.001.
Пример 6. Числа, представимые в ЭВМ, лежат в диапазоне При выходе результата за minD => underflow (исчезновение порядка), при выходе за maxD => overflow (переполнение).
•При переполнении обычно говорят, что плохи исходные данные, а при исчезновении порядка полагают результат равным нулю.
•Не следует торопиться. Пусть 10-78 |х| 1076 и вычисляется величина x = ab/(cd) при a = 10-30, b = 10-60, c = 10-40, d = l0-50.
•Если x=a•b/c/d => underflow; если x = 1/c/d•a•b => overflow.
•Если x = a/c•b/d => правильный ответ х = 1.
•Этот же ответ можно получить, если отмасштабировать переменные, например, умножив на 1040.
При переполнении или исчезновении порядка следует попытаться изменить последовательность действий, ввести масштабные множители и т. д.
При исчезновении порядка не всегда следует обнулять Пример 7. Пусть м = 10-2 и требуется найти S = 100 + 0.1 + … + 0. •Если вести суммирование слева направо, то с учетом округления до двух значащих цифр => S = •Если вычислять справа налево, то после тысячи слагаемых => 100, и дальнейшее {+0.1+0.1+…} ничего не изменит. Результат S = •Правильный результат S = 300 !? Как получить?!
•Сложим 1000 чисел по 0.1, затем еще 1000 чисел по 0.1, а потом сложим промежуточные суммы.
При сложении следует располагать слагаемые в порядке возрастания абсолютных величин, стараясь, чтобы при каждом сложении порядки величин различались мало.
При необходимости цикл суммирования разбивается на Аналогичное правило действует при перемножении большого числа сомножителей.
Последовательные приближения Пример 8. При расчетах методами последовательных приближений часто ведут вычисления до тех пор, пока поправка (разность между текущими и последующими приближениями) не станет меньше заданного порога. При этом, как правило, не обеспечивается заданная погрешность результата.
Если вести вычисления до тех пор, пока общий член ряда 1/k2 не станет меньше 10-3, т. е. до kобрыв = 32, и S = 1.610.
Правильно => S = 2/6 = 1,650....
Если же приближенную сумму ряда то погрешность останется бесконечной, как бы мала ни становилась Нужно помнить, что остановка итерационного процесса x1, х2, …. по косвенному критерию (например, по в задаче решения уравнения F (x) = O, не гарантирует достижения заданной погрешности Пример 9. Неустойчивость алгоритмов • Проверить неустойчивость алгоритмов (погрешность действия) на примере вычисления интеграла при помощи рекуррентной формулы • E0 вычислить аналитически и построить таблицу значений En при n = 1, 2,..., 20. Оценить возникающую ошибку.
• Повторить вычисления, изменив алгоритм на устойчивый En-1 = (1 – En) /n.
Аналитически и численно оценить ошибку при вычислениях En для n = 19, 18,..., 1 при выборе начального значения E20 = (показать, что начальная ошибка E20 < 1/21 и далее уменьшается !).
Пользуйтесь только устойчивыми численными Ю.Н. Прошин ЧМММ. Лекция 7 #
«