Кирьянов Д. - MathCad 11 (1077323), страница 50
Текст из файла (страница 50)
Для того чтобы ввести условие отражения на правой границе, пришлось определить еще один неизвестный пристрелочный параметр 220. Строка листинга, в которой определена функция score, задает условие стрельбы — сшивку двух решений вточке xf. В самой последней строке листинга выдан ответ — определенныечисленным методом значения обоих пристрелочных параметров, которыеобъединены в вектор и (мы применили в предпоследней строке операциютранспонирования, чтобы результат получился в форме вектора, а не матрицы-строки).Для корректного построения фафика решения лучше составить его из двухчастей — решения задачи Коши на интервале (xo,xf) и другой задачи Коши для интервала (xf,xi). Реализация этого способа приведена в листинге 12.4, который является продолжением листинга 12.3.
В последней строкелистинга 12.4 выведено значение второй искомой функции на правой границе интервала. Всегда полезно проконтролировать, что оно совпадает ссоответствующим пристрелочным параметром (выведенном в последнейстроке листинга 12.3).308Часть HI. Численные методыi Листинг 12.4.
Решение краевой задачи (продолжение листинга 12;3)V!м := 1010SO := rkfixedHoSI := rkfixedS0 M ) 2, 0 , 0.5 , M , D, 0.5 , 1 , M, DS1M.2 = 13-801Решение краевой задачи приведено на рис. 12.5. С физической точки зрения естественно, что интенсивность света уменьшается быстрее по мерераспространения в более плотной среде в правой половине расчетного интервала. В средней точке xf=o.5, как и ожидалось, производные обоих решений имеют разрыв.ПримечаниеЕще один пример решения краевых задач с разрывными коэффициентами ОДУприведен в справочной системе Mathcad.Рис. 12.5. Решение краевой задачи с разрывом в x f = 0 . 5(листинги 12.3—12.4)Ради справедливости необходимо заметить, что разобранную краевую задачулегко решить и с помощью функции sbval, заменив в листинге 12.2 зависимость а(х) на третью строку листинга 12.3.
В этом случае листинг 12.2даст в точности тот же ответ, что показам на рис. 12.5. Однако в определенных случаях (в том числе из соображений быстроты расчетов) удобнее использовать функцию bvaifit, т. е. вести пристрелку с обеих границ интервата.Глава 12. Краевые задачи309СоветЕсли Вы имеете дело с подобными уравнениями, попробуйте сначала решитьих как обычную краевую задачу с помощью более надежной и легкой в применении функции s b v a l .12.2. Задачи на собственныезначения для ОДУЗадачи на собственные значения — это краевые задачи для системы ОДУ,в которой правые части зависят от одного или нескольких параметров X.Значения этих параметров неизвестны, а решение краевой задачи существует только при определенных Х ь которые называются собственными значениями (eigenvalues) задачи.
Решения, соответствующие этим X*, называютсобственными функциями (eigenfunctions) задачи. Правильная постановка таких задач требует формулировки количества граничных условий, равногосумме числа уравнений и числа собственных значений. Физическими примерами задач на собственные значения являются, например, уравнение колебаний струны, уравнение Шредингера в квантовой механике, уравненияволн в резонаторах и многие другие.С вычислительной точки зрения, задачи на собственные значения оченьпохожи на рассмотренные выше краевые задачи. В частности, для многихиз них также применим метод стрельбы (см.
разд. 12.1.2). Отличие заключается в пристрелке не только по недостающим левым граничным условиям,но еще и по искомым собственным значениям. В Mathcad для решения задач на собственные значения используются те же функции sbval иbvaifit. В их первый аргумент, т. е. вектор, присваивающий начальныезначения недостающим начальным условиям, следует включить и начальноеприближение для собственного значения.Рассмотрим методику решения на конкретном примере определения собственных упругих колебаний струны.
Профиль колебаний струны у(х) описывается линейным дифференциальным уравнением второго порядка[pWdx I+ jLq^ytxJO,dx Jгде р(х) и q(x) — жесткость и плотность, которые, вообще говоря, могутменяться вдоль струны. Если струна закреплена на обоих концах, то граничные условия задаются в виде y(O)=y(i)=o. Сформулированная задачаявляется частным случаем задачи Штурма-Лиувилля. Поскольку решаетсясистема двух ОДУ, содержащая одно собственное значение X, то по идеезадача требует задания трех (2+1) условий.
Однако, как легко убедиться,уравнение колебаний струны — линейное и однородное, поэтому в любомслучае решение у(х) будет определено с точностью до множителя. Это оз-Часть III. Численные методы310начает, что производную решения можно задать произвольно, напримеру (0)=i, что и будет третьим условием. Тогда краевую задачу можно решатькак задачу Коши, а пристрелку вести только по одному параметру — собственному значению.Процедура поиска первого собственного значения представлена в листинге 12.5.Листинг 12.5. Решение задачи о собственных колебаниях струныр (х):= 1a := Оа (х):= 1р ' (х):= ОЬ := 1Хо := 0 .
5У1D ( х , у):=•(р'( х ) - y i + у г ' q ( x ) •;Р (х;(0l o a d (a , Я) :-s c o r e (Ь , у) := уоЛ :- s b v a l \к , a , b , D , l o a d , s c o r e )Л - ( 9.871 2 - я 2 = 9.87В первых двух строках листинга определяются функции, входящие в задачу,в том числе р ' ( х ) : = о , и границы расчетного интервала (0,1). В третьейстроке дается начальное приближение к собственному значению А„, в четвертой вводится система ОДУ.
Обратите внимание, что она состоит не издвух, а из трех уравнений. Первые два из них определяют эквивалентную (I)систему ОДУ первого порядка, а третье необходимо для задания собственного значения в виде еще одного компонента у2 искомого вектора у. Поскольку, по определению, собственное значение постоянно при всех х, тоего производная должна быть приравнена нулю, что отражено в последнемуравнении. Важно также, что во втором из уравнений собственное значениезаписано как у2, поскольку является одним из неизвестных.В следующих двух строках листинга задается левое граничное условие,включающее и недостающее условие на собственное значение для третьегоуравнения, и правое граничное условие у0=о.
В предпоследней строке листинга обычным образом применяется функция sbval, а в последней выводится результат ее работы вместе с известным аналитически собственным2 2значением п я . Как легко убедиться, мы нашли первое собственное значе-Глава 12. Краевые задачи311ние для n=i, а чтобы найти другие собственные значения, необходимо задать другие начальные приближения к ним (в третьей строке листинга 12.5).Например, выбор Х0^50 приводит ко второму собственному значению 2 2 л 2 ,а >.0=юо — к третьему з2-п2.Чтобы построить график соответствующей собственной функции, надо добавить в листинг строку, программирующую решение задачи Коши, например, такую: u : =rkfixed(ioad(a,A) ,a,t>,ioo,D).
Полученные кривые показаны на рис. 12.6 в виде коллажа трех графиков, рассчитанных для трехсобственных значений.Примечание)Примеры решения нескольких задач на собственные значения можно найти вразделе Mathcad Resources.Рис. 12.6. Первые три собственные функции задачи колебаний струны(коллаж трех графиков)12.3. Разностные схемы для ОДУМногие краевые задачи не поддаются решению методом стрельбы.
Однаков Mathcad 11 других встроенных алгоритмов нет. Тем не менее, это не означает, что по-другому решать краевые задачи невозможно, ведь другие численные алгоритмы несложно запрограммировать самому пользователю. Рассмотрим возможную реализацию наглядного метода, называемогоразностным, которым можно решать краевые задачи как для ОДУ, так и длядифференциальных уравнений в частных производных.12.3.1. О разностном методе решения ОДУРазберем идею разностного метода решения краевых задач на примеревзаимодействия световых пучков (см. рис. 12.1), переобозначив в системе(12.1,1) интенсивность излучения вправо на у, а интенсивность излучения11 Зак 984312Часть III.
Численные методывлево на у (просто в целях удобства, чтобы не писать индекс). Суть методазаключается в покрытии расчетного интервала сеткой из ы точек. Тем самым определяются ( N - D шагов (рис. 12.7). Затем надо заменить дифференциальные уравнения исходной краевой задачи аппроксимирующими ихуравнениями в конечных разностях, выписав соответствующие разностныеуравнения для каждого i-ro шага. В нашем случае достаточно просто заменить первые производные из (12.1,1) их разностными аналогами (такой метод называется еще методом Эйлера):Д(1)-=-i-=^ f v . l • \/•— W v .