Форсайт Дж., Малькольм М., Моулер К. - Машинные методы математических вычислений (1040536), страница 6
Текст из файла (страница 6)
Предположим, что лишь для одного нз двадцати коэффициентов делается ошибка в младшем двоичном разряде. Именно, предположим, что коэффициент при х" изменяется Е„= О.О, Е„0.0500000, Е~о . 00500000 Е„0.0527778, Е„0.0557190, Ем 0 0590! 76, Е„0.0627322, Еь о = 0 0669477, Еоо 0 0717733, Е,о 0 0773523, ЕРо 0.0838771, Ео 0.0916123. а Вычисления с пллВлюшеи тОчкОЙ 8 32 10,09526 6145 х- 0,64350 09041 11.79363 388! с- !.65232 9728! 13.99235 8137 ~ 2.51883 0070с 16.73073 7466 ~ 2.81262 4894с 19,50243 9400 ~ 1.94033 0347! 1.00000 0000 2.00000 0000 3.00000 0000 4.00000 0000 4.99999 9928 6.00000 6944 Г>.99969 7234 8.00726 7603 8.91725 0249 20Д4690 8!01 Заметим, что малое изменение в коэффициенте †2 имело следствием то, что десять корней стали комплексными и что два из ннх отодвинулись от действительной оси более чем на 2.81.
Конечно, полная запись р(х) в машине потребовала бы значительно болыпе округлений, а вычисление корней сопряжено с еще ббльшими ошибками. Приведенная выше таблица корней была результатом очень точных вычислений и не пострадала сколько-нибудь заметно от ошибок округления.
Причина, по которой эти корни так сильно изменились, лежит ие в ошибках округления и не связана с выбором алгоритма вычисления; суть в чувствительности самой задачи. Легко провести анализ того, что произошло. Можно записать поливом в виде р (х, а) = — х" — ах" +... и затем найти частную производную по а для каждого корня полинома р(х). Зто делается дифференцированием уравнения р(х, а)=0 по а: др (х, сх) дх др (х, а) дх да да дх дрсда хсо дР~дх оо оо Ц (х — () с=с с= ! / Вычисляя зто выражение для каждого корня, получаем дх! см — — с = 1, 2, ..., 20.
дсх ~„; оо П (' — ') / =с (~с с — 210 на — 210+2 . Какое воздеиствие зто малое пзменениес произведет иа корни полинома? Чтобы ответить на этот вопрос, Уилкинсон тщательно вычислил (при р=2, 1=-90) корни уравнения р(х)+2 "х"=-О. Приводим их значения, правильно округленные до указанного числа цифр: ЕО. РЕШЕНИЕ КВАДРАТНЫХ УРАВНЕНИЙ Эти числа дают прямую меру чувствительности каждого нз корней к изменению коэффициента а. Приводим их значения: ах~а«1.. Корень 2.6. Ре1иеиие ввадйоатнык уравнензей Имеется знаменитый алгоритм для решения квадратных уравнений, заключенный в следующей математической теореме: Теорема. Если а, Ь, с — действительные числа и аФО, то уравнение ах*+ох+с=О имеет ровно два решения, именно — Ь+ 1/Ь вЂ” 4ве 2а — Ь вЂ” Р' Ь» — 4оо х» = Посмотрим, как работают эти формулы, будучи применены «в лоб» для вычисления х, и х,, На этот раз мы возьмем плавающую систему с 15=10, 1=8, Е= — 50, й'=50; она точнее многих широко используемых машинных систем.
1 2 з 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 2о 8 2 х 10-~« 8.2 х 1О- -1.6 х Ш-« 2,2х10» -6.1 х 1О « 5,8 х 10» -2.5 х 10» 6.О х ГО -8.1 х 1О» 7.6 х 1О« —.4.6 х 10« 2.0 х 10» 62 х 10« 1.З х 1О -2Л х 10' 2.4 х 10» -1.9 х 1О 1.О х 1О -Зй х 10« 4.1 х 1О 34 2. ВЫЧИСЛЕНИЯ С ПЛАВАЮЩЕЙ ТОЧКОЙ Случай 1. а=1, Ь= — 10', с==1. Точные корни соответствующего квадратного уравнения, правильно округленные до 11 значащих десятичных цифр, суть х, = 99999.999990, х, = 0.000010000000001. Если мы воспользуемся формулами теоремы, то вычислим х, = 100000.00 (очень хорошо), х,=-0 (совершенно неверно). (Советуем читателю проверить, понимает ли он, как х, стано- вится нулем в этом плавающем вычислении.) Снова, на этот раз при вычислении х, мы пострадали от катастрофической потери верных знаков.
Имеются различные альтернативные способы вычисления корней квадратичного полинома, которые избегают подобной потери. Одним из такрх способов является использование знака д для определения, какая формула влечет меньшую потерю знаков, и затем вычисле- ние одного из корней по этой формуле. Таким образом, Ь+ Меп (Ь) Г' Ь~ — 4ас х =-— 1 2а Поскольку из равенства ах'+эх+с=а(х — х,) (х — х,) следует ах,х,=с, то второй корень можно вычислить по формуле с х = —.
ек, ' Для данных случая 1 такой способ дает х,=100000.00, х,=1.00000001100000.00=-0.000010000000; оба результата вполне приемлемы. Сейчас мы хотели бы предложить следующий критерий для реализации машинного алгоритма решения квадратных уравнений. Он формулируется здесь довольно вольно, однако в статье Форсайта (!969) можно найти строгую формулироику Мы скажем, что комплексное число г находится строго в пределах множества Р, если либо г=0, либо Ис + 4 < Яе г =- (1и - г и ()с" < 1ш г С )) и - з.
Это означает, что действительная и мнимая части г находятся строго внутри области чисел, которые могут быть хорошо приближены элементами из Е, Довольно произвольный множитель р' привлечен для того, чтобы обеспечить некоторый запас надежности. Предположим, что все три числа а, Ь, с находятся строго в пределах Р. Тогда они должны быть допустимы в качестве Ез. РЕШЕНИЕ КВАДРАТНЫХ УРАВНЕНИЙ зз входных данных для алгоритма решения квадратных уравнений.
Если а=-Ь= — с=О, то алгоритм должен заканчивать свою работу сообщением, что все комплексные числа удовлетворяют уравнению ах'+Ьх+с=О. Если а=-Ь=О и с~О, то на выходе алгоритма должна быть информация о том, что ни одно комплексное число не удовлетворяет этому уравнению. В противном случае пусть г, и гз — точные корни уравнения, занумерованные так, что !гз~(1гзй (Если а=-О, положим гз=оз.) Во всех случаях, когда г, находится строго в пределах множества г", алгоритм должен находить хорошее приближение к еб под этим понимается требование, чтобы приближение отличалось от г, не более чем, скажем, на одну единицу предпоследнего знака корня.
То же самое должно выполняться для ез. Если один или два корня не находятся строго в пределах множества г, то должно выдаваться соответствующее сообщение и тот корень, который лежит строго в пределах г (если таковой имеется), должен быть определен с указанной точностью. Зтим завершается вольное описание желаемого исполнения алгоритма для решения квадратных уравнений. Вернемся теперь к рассмотрению некоторых типичных уравнений, чтобы посмотреть, как работают для них квадратичные формулы. Случай 2. а=б, Ь=5, с= — 4.
Здесь нет някаких затруднений при вычислении х,= =0.50000000 и хз=- — 1.3333333 или приближений к этим значениям, какой бы формулой ни пользоваться. Случай 3. а — -бх!0*', Ь=5х!О", с= — — 4х10зз. Поскольку, с точностью до множителя 1О", это те же коэффициенты, что н в случае 2, то корни не изменились. Однако, применение формул для х, и х, вызовет переполнение, так как Ь', будучи больше 10", находится за пределами г.
Скорей всего, этн равномерно большие значения для 1а1, !Ь1, 1с1 смогут быть обнаружены до входа в алгоритм, и все три числа будут поделены на какой-нибудь масштабирующий множитель типа 10", после чего задача сведется к случаю 2. Случай 4. а=10-зз Ь вЂ” 10зз с=!Озз Здесь х, примерно равен 1, в то время как х, приблизительно равен !О". Следовательно, наш алгоритм должен определить х, очень точно, несмотря на то что х, находится за пределами Е. Очевидно, что всякая попытка достигнуть примерного равенства модулей коэффициентов простым делением их иа одно и то же число обречена па неудачу и сама по себе может повлечь переполнение илн появление машинного нуля.
Зто уравнение в самом в Вычисления с плАВАющея точкОЙ деле является суровым испытанием для любой программы, решающей квадратные уравнения. Читатель может подумать, что квадратное уравнение, у которого один корень находится внутри пределов Р, а другой — ! вне их, есть искусственный пример, не имеющий практического; значения. В таком случае он ошибается. Во многих итерационных алгоритмах, содержащих подпрограмму решения квадрат-. ных уравнений, квадратичные многочлены имеют вырожденное поведение, характеризуемое тем, что а- 0 в случае сходимостн алгоритма. Пример такой ситуации — метод, предложенный в статье Мюллера (1956), для нахождения нулей гладкой функции общего вида.
Случай б. а=!.0000000, Ь= — 4.0000000, с=-3.9999999. Здесь корнями являются х,=-1.999683?72, х,=2.000316228. Однако применение квадратичных формул дает х, =-х,=2.0000000, и только первые четыре знака верны. Эти приближения далеко не удовлетворяют выдвинутому выше критерию, однако по сравнению с другими примерами трудность здесь совсем иная. Данный квадратичный полинам соответствует квадратному уравнению (х — 2)'=е, где е=0.000000!. !(ак уже говорялось в предыдущем параграфе, это случай чувствительности самого уравнения, а не метода его решения.
Однако чувствительные квадратные уравнения все же могут решаться так точно, чтобы удовлетворить нашему критерию, если часть вычислений проводить с повьппенной точностью (в данном случае должно быть 1=16 или еще больше). В определенном смысле корни, вычисленные в случае 5(х,=х,=2), являются «хорошими» решениямн этого уравнения. Заметим, что они будут точными корнями уравнения х' — 4х+4 = О. Это то же уравнение, что и в случае 5, за исключением свободного члена, отличающегося от с единицей последнего разряда, Таким образом, мы вычислили точные корни уравнения, очень «близкого» к тому, которое пытались решить.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
