Говорухин В., Цибулин Б. Компьютер в математическом исследовании (1185927), страница 76
Текст из файла (страница 76)
15.14: » сапет(т(:. 1).т(:,3)) Ту же фазовую траекторию можно нарисовать при помощи команды » 01от(т(:.1)Л(:.3).'И'.У(еп0,1).т(епп,з).'ои') Решение краевых задач 407 Различные примеры приготовленных правых частей содержатся в файлах огЬ)гоое (ограниченная проблема трех тел), огЬс2осГе (задача двух тел), пдЫосГе (уравнения Эйлера), чг)росГе (уравнение Ван дер Поля). 40 ЗО то О -20 -то о то 20 Рис. 1В.14. Фазовая травкторня, выведенная командой соввт В заключение приведем еще одну форму обращения к команде осГе 5ОГ нег решения дифференциального уравнения: ГТ.У.5)-все во1тег("Р ,Т5РАИ,ТО,ОРТГОЙ5,РАРАИ5) По сравнению с ранее рассмотренными обращениями здесь в левой части появился параметр 5, в котором накапливается статистика о ходе решения задачи. Зто массив из шести элементов, в которых содержится следующая информация: о число успешных шагов; О число неудачных попыток; о число вычислений правой части; о число вычислений матрицы Якоби; о число проведенных факторизаций; сз число решений систем линейных алгебраических уравнений.
решение краевых задач В МАТЮКАВ 6 добавлена команда Ьчр4с для решения одномерных краевых задач, задаваемых системой обыкновенных дифференциальных уравнений первого порядка. Решение краевой задачи на отрезке ищется при помощи метода коллокаций, для дискретизации используются формулы Лобатто и конечные разности. Как обычно в МАТВАВ, действуют разумные назначения точности, допустимого числа разбиений и ряда других параметров.
Для получения их значений имеется команда Овраг, а изменить их позволяет функция Ьчрвег. Кроме того, имеются вспомогательные функции для формирования начального приближения (Ькрт От 1) и получения решения в заданных точках (ЬчркзГ ). Для использования функции Ьчр4с краевая задача должна быть представлена в виде системы обыкновенных дифференциальных уравнений первого порядка относительно неизвестной вектор-функции у(х) переменной х и параметров р: 408 Глава 15.
Численный анализ в МАТЕАВ у'-((х,у,р) Задача рассматривается иа отрезке [а. Ь3, и краевые условия формулируются следующим образом: с(у(а) . у(Ы, р)-0 Здесь с — заданная функция от значений функций в точках а и Ь, Обращение к Ьчр4с имеет следующий вид: 50Е=Ьчрлс(ООЕ, ВС. 1М 11. ОРТ. Р1. Р2...
) Здесь ООЕ и ВС обозначают имена функций, предназначенных для вычисления правой части системы и краевых условий соответственно. Структура 1М1Т с зарезервированными полями х и у дает начальное приближение к решению, узлы идут в порядке возРастания, так что 1М!Т. х(1)=а и 1М1Т. х(еп0)=Ь. Структура 1МП может быть сформирована при помощи команды Ьчр)0)Ь.' Ниже в примере указаны два способа формирования структуры, задающей начальное приближение. Имя ОРТ обозначает переменную с параметрами расчета, а параметры задачи Р1, Р2, ... должны передаваться в функцию вычисления правой части и функцию, задающую краевые условия.
Способ передачи аналогичен рассмотренному ранее для команд решения задачи с начальными данными. Если в системе имеются параметры Р1, Р2,..., то место параметров ОРТ при обращении к Ьчр4с должно быть занято хотя бы пустым массивом П. Если интегрируется система без параметров и иет необходимости менять параметры расчета, то обращеиие к функции Ьчр4с может выглядеть следующим образом: 50( Ьчрас(ООЕ.ВС.!М1Т) Рассмотрим оформление заголовков функций вычисления правой части ООЕ и задания краевых условий ВС. Входными параметрами для функции вычисления правой части с именем ООЕ являются точка отрезка Х, вектор значений в этой точке Т и параметры задачи Р1, Р2,.... Выходные данные представлены вектором ОТ: ОТ ООЕ(Х,Т,Р!,Р2,...) Краевые условия обрабатываются функцией С-ВС(та,ТЬ,Р1.Р2....) Здесь векторы Уа и УЬ дают значения вектора иа разных концах интервала, а выходиым параметром будет вектор С.
Решение 50Е вычисляется в виде структуры, поля которой фиксированы и имеют следующий смысл; х — координаты узлов разностиой сетки, у — аппроксимация решения в узлах, ур — аппроксимация первой производной решения в узлах, и рагаветегз — значения неизвестных параметров. Последний параметр позволяет решать задачи по определению собственных чисел и функций (см, файл л)а14Ьчр. л) из каталога демонстраций фю(Ьох~л)ас(аЬ~делзоз, где Ьчрйс применяется для вычислеиия четвертой собственной функции уравнения Матье). Для изменения характеристик расчета можно задать соответствующие параметры ОРТ: АЬВТОТ дает величину абсолютной, а Ве1То! — относительной погрешности, параметр ЕОасоЬ) ал указывает иа наличие аналитического выражения для якобиаиа системы, а ВСОасоЬ)ал — якобиаиа для системы краевых условий, Мяах задает предельиый порядок получающейся системы линейных алгебраических уравнений, Решение краевых задач 409 50асз определяет параметр решения.
Порядок системы есть число неизвестных функций, умноженное на число узлов, по умолчанию число )йвах равно 1000. Ко- манда Ьчрзег() распечатывает имена всех параметров и их возможные значения, отмечая используемые по умолчанию фигурными скобками. В процессе расчета происходит измельчение сетки для достижения заданной точности. Рассмотрим пример вычисления решений следующей краевой задачи: уееуг 0 у(0) 0 у(2) .1 Оформим пример с использованием подфункций зоьтипсс! оп.
В одном файле поме- стим операторы подготовки н вычисления решения, а также вспомогательные функ- ции вычисления правой части обе, задания краевых условий Ьс и подготовки на- чального приближения топ. При передаче в качестве параметров имен подфункцнй снабжаем их префиксом Ф. Так, при первом вызове функции Ьчр1 п)0 испо)) рзуем функ- цию топ для задания начального приближения решения в виде синусойды и произ- водной решения в виде косинусоиды.
Для построения второго решения задаем на- чальное приближение константами для каждой компоненты решения. Затем дважды вызываем функцию решения краевой задачи Ьчр4с для нахождения решений н ис- пользуем команду Ьчрча) для получения массивов решений в заданных точках. Гопсг!оп Ьчр бекю Хтно зо)иг!оп Ь"2; п 8; х Ып*[0:п); Чшзг-Ьчрзп11(х.ФПзп): Ч!=Ьчр4с(Фобе,ФЬс.у!п!С): Ч!пз(-Ьчр1шг(х,[-1 01); Ч2-очр4с(Фобе.ФЬс,у!п!1); х Ы50*[0:501; у1-Ьчрча)(Ч1.х); у2-Ьчрча)(Ч2.х): зчор)ог(121).
р)ог(х.у1(1.:),х.у2(1.:)); зчЬр)о((122) . р) ог(х,у!(2.: ), х.у2(2,: ) ); (ипсо!оп Ч-Гоп(х) Ч"4*[51п(2.*х): соз(3.*х)): топот!оп бубх-обе(х.у) бубх-[у(2): -у(1)."21: гипс(!оп Г-Ьс(уа.ую 1-[уа(1): УЫ1)+Ц; ъ,о 15 2 О 0.6 -1 О 0.5 1 х 1 1.6 2 х Рис. За.зй. Дее решение краевой задачи 410 Глава ть. Численный анализ в ИАТЕАВ Результат выполнения функции Ьзр сизо представлен на рис.
15.15, в левой части даны сами решения — функции у(х), а на правом рисунке приведены графики их производных. Отметим, что для вычисления первого решения (пунктир) потребовалось 19 узлов сетки, включая концы интервала, а для вычисления второго решения (сплошная кривая) хватило 10 узлов. Решение начально-краевых задач параболического типа В МАТ).АВ 6 появилась также команда роере для решения системы уравнений в частных производных параболического типа, с одной пространственной координатой. Общий вид системы относительно неизвестной вектор-функции Ет[х О может быть представлен следующим уравнением: С вЂ” = х — (л'"г )+ 6 а~ Вх Здесь матрица С, векторы г и 3 зависят от переменных х СЕТ В .
Матрица С является диагональной, по крайней мере, один элемент ее должен быть ненулевым. Для уравнения теплопроводности вектор Г определяет поток, а вектор 6 описывает источники. Введение параметра лз позволяет решать одномерные задачи для радиальной координаты в цилиндрической или полярной (т-1), а также сферической (и 2) системах координат. Начальные условия определяются для вектора неизвестных: ЕГ(х,тв) = Е)л(х) Краевые условия задаются в следующем виде: Р(хдИ+ Цх ОГг(х ге(,— ) = О аи 'ах Обращение к функции решения системы уравнений в частных производных имеет вид: 50$.-рнера(в. РОЕ,1С.ВЕ, ХХ,ТТ.ОРТ, Р1, Р2, .. Л Здесь параметр в определяет систему координат: О отвечает декартовым координатам, 1 — цилиндрическим, 2 — сферическим.
Вычисление правой части рассматриваемой системы производится в функции РЮЕ, данные о начальных условиях берутся из функции 1С, а краевые условия задаются функцией ВС. Сетка по координате х на интервале [а, Ь) должна содержать не менее трех узлов и задается массивом ХХ, содержащим узлы в порядке возрастания координаты. Лучше использовать более мелкий шаг между узлами в тех областях, где решение меняется значительно, и более крупный шаг в местах плавного изменения решения. Для в О (цилиндрические или сферические координаты) не обязательно сгузцать сетку в окрестности особой точки х-О, и при в О значение а не может быть отрицательным, Для интегрирования задачи по времени в функции рс)ере используются команды для решения задачи Коши (опе1бв и др.), а применяемый метод и шаг ннтегрирова- Решение начально-краевых задач параболического типа 411 ния выбираются функцией рбере.
Массивом ТТ определяется набор временных слоев, в которых будут запомнены решения, таких слоев должно быть не менее трех. Некоторые расчетные параметры можно изменить при помощи команды обезеС, см. описание в разделе, посвяшенном интегрированию задачи Коши. В частности, доступны такие параметры, как относительная (Ое1ТЬ1) и абсолютная (АЬСТо)) величины погрешности, начальный (!п)С1а15Сер) и максимальный (йах5Сер) шаги интегрирования.