Говорухин В., Цибулин Б. Компьютер в математическом исследовании (1185927), страница 93
Текст из файла (страница 93)
1); бе(М)-Ь; бас(:,К)-(теча1(тцп.а+бе)-тг).т'Ь; епб епб Функцией ]асо завершается файл печ(тона.(п. Для решения конкретной системы нужно подготовить файл, в котором будет вычисляться вектор-функция и, если надо, матрица Якоби. В качестве примера п)-файла, задающего уравнения системы и матрицу Якоби, приведем функцию, позволяющую вычислять одновременно несколько корней Рецтение нелинейных уравнений 493. кубического полинома.
В этом случае каждая компонента вектор-функции зависит только от соответствующей компоненты вектора неизвестных и матрица Якоби диагональная: гипс(топ [г,зас.пац]-гоп(г,па9) тт лагун<2 > тяеаргу(Пац> ) -Па9. Г!ац-0; )ас"П: епо Г г. 3+1: >Г П ац, Зас=.отед(3 г."2); епо Представим результат вычисления сразу трех корней. Отметим, что такой способ возможен только для нахождения комплексных корней скалярного уравнения. Итерации производятся до тех пор, пока условие сходимости не выполннтся для всего вектора решения: » [х,п]=пеи(опа(>топ",[1 т -т].20.[],[],'геяц>(дхс') х -1.0000 0.5000 - 0.8660т 0.5000 + 0.8660т и- 10 Этой функцией можно пользоваться и при вычислении одного корня. Проверим заодно, что получается, если число итераций невелико: » [х,п]-пеитопя("с05 Г)3",7.5.[],[],'геац>С.Схг'): Нет схоаииости 5 итераций Теперь посмотрим результаты, помещенные в файл геен[с,схт: Итерация Нарна > Решение 1 1.41Е+000 О.ОООЕ»000»т(-1.0ООЕяООО> 2 5.97Е-001 З.ЗЗЗЕ-001+т(-6.667Е-001) 3 3.31Е-О01 5.822Е-001+т(-9.244Е-001) 4 2.73Е-002 5.088Е-001+7(-8.682Е-001) 5 2.44Е-004 5.001Е-001+1(-8.660Е-001) б 1.98Е-О08 5.0ООЕ-001т-т(-8.660Е-001) Нет сходииости 5 итераций Видно, что для достижения точности и определения корня не хватило всего одной итерации.
В общем случае матрица Якоби не обязательно диагональная. Например, если записать кубический полипом в виде двух уравнений относительно вещественной и мнимой частей неизвестной, то и>-файл будет выглядеть следующим образом: Гцпсгтоп [Г,]ас. П ац]-[сп(х, Пац> т' - [х(1)*(х(1)"2-Зих(2)12)+1. х(2)*(3*х(1)"2-х(2) 2>]; тт Г1ац Зас [3*х( 1) 2-3"х(2) 2 -6»хы )*х(2); бих( 1)*х(2) 3*х( 1) 2-3'х(2) 2]: епц Если матрицу Якоби находить численно, то функция для определения корней кубического полинома будет иметь простейший вид: гцпсстоп Г-Гоп(г) Г а',»841; 494 Глава 1в. Дополнения и примеры Бассейны для корней кубического полинома Рассмотрим кубическое уравнение 2м1-0, имеющее один действительный и пару комплексно-сопряженных корней.
В зависимости от начальной точки метода Ньютона будут сходиться к разным корням, Область начальных значений, из которых итерации метода Ньютона приводят к одному и тому же корню уравнения Г(х) О, назовем бассейном корня, Для рассматриваемой задачи комплексная плоскость делится на три части по числу корней, и при этом границы бассейнов имеют фрактальную структуру. Подготовим программу рисования бассейнов для корней кубического уравнения. Будем окрашивать бассейн каждого корня своим цветом, а для большей изобразительности используем изменение интенсивности цвета в зависимости от числа итераций, которые пришлось вьшолнпть для нахождения корня с заданной точностью. Быстрой сходимости к решению будет отвечать темный тон, и интенсивность будет убывать по мере увеличения числа итераций.
Подготовим данные для расчета бассейнов в области 1-3,1]х(-1,5,1.51 (параметры а1,а2,01,Ь2) с числом точек 400х300. Введем массив начальных точек 70, массив результатов 77 н массив для хранения числа потребовавшихся итераций ЬЬ: » а1--3; а2=1. Ь1=-1.5: Ь2=-Ь). и-200; и-150. Пх-(а2-а1) Г(2'п-1>; Пу-(Ь2-Ь1) Г(2*». 1): (Х.У]-аеапог1а(а1:Ьх:а2.Ы:Ьу Ь2): 70-Х;*У; е1- 001; ГЬ-еегоа(2"е,2"и); 77-Е, Запустим вьпшслеиия методом Ньютона для всех точек из массива 70. Для данного примера не будем пользоваться обращением к приведенной выше функции пеыбопв, а запишем явиуку формулу для вычисления итераций: » Гог и=1:)еп9гп(70(:>).
ХО-70(Х); а-ХО»1: ию)е ава(а-гО»е1. ЬЬ(Х) ШХ)+1: 40-г; а-(2»а"5-1)Г(2»г"2); епа, 27(Х)-г; епа Для подготовки рисунка преобразуем данные, смягчив логарифмированием отличия по числу потребовавшихся итераций и масштабировав элементы массива ЬЬ от нуля до единицы: » Ь~.-)о0()оО(ЬЬ)~1); Ьап и в)п(и>п(ЬЬ)): Ьиах.иах(иах(ЬЬ)); ЬЬ ((с-Ьи1п)Г((иах-!и1п); Подготовим палитру из 66 интенсивностей серого, удалив темные и близкие к белому тона: » пс-бб; сс-опеа(пс,З); Гог Х-1:пс; сс(Х,:)-.45 ХГпс*,45: епа; со) огиар(сс); Вычислим массив С, используя для каждого корня свой диапазон интенсивности серого тона.
Обратим внимание, насколько просто в МАТЬАВ произвести разделение начальных точек по корням, к которым привели итерации метода Ньютона. Решение нелинейных уравнений 495 1.5 0.5 -0.5 -1.5 '-з Рис. 1ад. Бассейны корней. Серая палитра с монотонным убыванием интенсивности 0.5 -0.5 -1.5 -3 -1 0 Рис. 1а.Я. Серая палитра с вариациями интенсивности Затем обратимся к команде )рщде для вывода изображения, подготовленного массивом С. На рис. 18.3 приведен результат.
» С-т1Х(1~(СС+(таад(22)>Е1) -2и(1аад(22) -Е)) )ЯПС/3), 1ма де( Х(1,, ) . У(; . 1) . С ), а х1 в (ваде, акт в( < а! а2 Ь) Ь2) ) Усоверщенствуем рисунок, изменив палитру так, чтобы интенсивность каждого второго тона чуть возросла (напомним, что нуль отвечает черному цвету, а единица — белому). Для этого достаточно выполнить следующую команду: » топ и-2;2:пс. сс(Н.:)-сс(Н.;)*.95;епО, со)опмр(сс); Полученное в результате изображение приведено на рис. 18.4.
496 Глава 18. Дополнения и примеры Изменим параметры области так, чтобы рассмотреть распределение точек в окрест- ности нуля: » а1--1: а2-.5: Ы--.25; 02--01: и-240; йт.60; Полученное изображение приведено на рис. 18.5. Хорошо видно, что множество начальных точек, из которых метод Ньютона сходится к данному корню, имеет фрактальный характер. 0.25 -0 25 -1 -0.5 0.5 Рик. 18.5.
Фрагмент картины бассейнов корней Для сравнения приведем результаты вычисления бассейнов при помощи процедур пакета ХАС. В левой части рис. 18.6 дана карта бассейна для корней кубического полинома, рассчитанных при помощи функции с05пЬГ, а в правой части— соответственно для с05ОЬГ.
Каждому корню отвечает своя интенсивность серого, а белым цветом даны точки, итерации из которых закончились сообщением об аварийном завершении работы процедур. Видно, что алгоритмы из библиотеки НАС дают картины, отличные от полученных прямым методом Ньютона, см. рис. 18.3. В отличиие от рис. 18.3 центром сгущения областей является начало координат.
Таким образом, результаты численного исследования бассейнов кубического полинома существенно зависят от используемых методов. -3 -2 -1 0 1 -3 -2 -1 О 1 Рис. 18.6. Бассейны корней для процедур из библиотеки МАБ Для получения цветного рисунка можно воспользоваться одной из многих палитр МАТ1.АВ, см.
главу 14 «Графика МАТВАВ», или подготовить собственную. Например, палитру из трех цветов (синий, зеленый и красный) с оттенками можно создать при помощи следующих команд: » пс-66; пс1-11; сс-кегов(пс.З>; Гог К-1;пс1, Разработка приложения с 6))1 497 сс(к.З)-.5 К/пс1*.5: сс(К+пс1,3)-1; сс(к+пс1.2)-.4+К/пс1*.6; сс(К+2 пс1.2)-.4+к/пс1*.3; сс(К»З*пс1,2)-0.7»К/пс1*.3; сс Р»4"пс1. 1)-. 5»К/пс1*.
4; сс(К+5*пс1,2]-.5 к/пс1*.4: сс(К+5*пс1.1)-1: епо' Разработка приложения с 6Ш Опишем последовательность действий, необходимую для создания интерактивного графического приложения. Представим пример еручного» программирования графического интерфейса и приведем тексты на языке МАТВАВ, реализующие соответствующие конструкции. Конечно, используя визуальное проектирование при помощи СпЫе, решать подобные задачи проще, но вместе с тем полезно знать, как именно устроены команды, реализующие поддержку графики, и опыт такого рода может помочь при отладке графических приложений н квалифицированной доводке сложных программ.