Кирьянов Д. - MathCad 11 (1077323), страница 55
Текст из файла (страница 55)
Таким образом, зависимость от найденного напредыдущем шаге по времени решения определяется зависимостью от негоисточника, т. е. правой части.Суммируя сказанное, можно констатировать, что если Вы имеете запрограммированный алгоритм решения выписанного эллиптического уравнения, чуть более сложного, чем уравнение Пуассона, то его с легкостьюможно использовать в качестве подпрограммы реализации разностной схемы двумерного уравнения теплопроводности. Забегая вперед, приходитсяотметить, что, к сожалению, встроенные функции Mathcad для решенияуравнения Пуассона в данном случае не годятся в качестве такой подпрограммы, поскольку предполагают независимость источника от самой неизвестной функции и могут справляться лишь с правой частью, зависящейтолько от пространственных координат (см.
разд. 13.3.2).Не будем далее останавливаться на способах решения многомерных уравнений, ограничившись этими замечаниями относительно путей оптимизацииалгоритмов их решения.13.3. Встроенные функции для решенияуравнений в частных производныхКак видно из предыдущего раздела, с уравнениями в частных производныхвполне можно справиться, и не прибегая к специфическим средствамMathcad. Между тем, в пакете Mathcad 11 имеется несколько встроенныхфункций, при помощи которых можно автоматизировать процесс решениядифференциальных уравнений в частных производных.
Рассмотрим в данном разделе основные аспекты их применения, отмечая не только инструкции по их применению, описанные разработчиками Mathcad, но также инекоторые «скрытые» возможности этих функций.13.3.1. Параболические и гиперболическиеуравненияВ новой версии Mathcad 11 разработчики впервые применили встроеннуюфункцию pdesoive для решения уравнений в частных производных, отлично осознавая значимость этих задач для современного исследователя и инженера. Эта функция применяется в рамках вычислительного блока, начинающегося ключевым словом Given и пригодна для решения различныхгиперболических и параболических уравнений.Встроенная функция для решения одномерного уравнения (или системыуравнений) в частных производных (того, которое определит пользователь в338Часть ill.
Численные методырамках вычислительного блока Given), зависящего от времени t и пространственной координаты х, имеет целый набор различных аргументов иработает следующим образом.ПPdesolve(u, x, xrange, t, trange,[xpts] ,[tpts] ) ) — возвращаетскалярную (для единственного исходного уравнения) или векторную(для системы уравнений) функцию двух аргументов ( x , t ) , являющуюсярешением дифференциального уравнения (или системы уравнений) вчастных производных. Результирующая функция получается интерполяцией сеточной функции, вычисляемой согласно разностной схеме.• и — явно заданный вектор имен функций (без указания имен аргументов), подлежащих вычислению.
Эти функции, а также граничныеусловия (в форме Дирихле или Неймана) должны быть определеныпользователем перед применением функции pdesoive в вычислительном блоке после ключевого слова Given. Если решается не системауравнений в частных производных, а единственное уравнение, то, соответственно, вектор и должен содержать только одно имя функции ивырождается в скаляр.• х —пространственная координата (имя аргумента неизвестной функции).• xrange — пространственный интервал, т. е.
вектор значений аргумента х для граничных условий. Этот вектор должен состоять из двух действительных чисел (представляющих левую и правую границурасчетного интервала).• t — время (имя аргумента неизвестной функции).•trange -— расчетная временная область: вектор значений аргумента t,который должен состоять из двух действительных чисел (представляющих левую и правую границу расчетного интервала по времени).•xpts — количество пространственных точек дискретизации (может неуказываться явно, в таком случае будет подобрано программой автоматически).•tpts — количество временных слоев, т. е.
интервалов дискретизациипо времени (также может не указываться пользователем явно).ПримечаниеПомимо этой функции для решения параболических и гиперболических уравнений, начиная с новой версии Mathcad 11, можно использовать еще одну встроенную функцию numol ( ) . Функция numol () имеет еще большее число аргументов и позволяет управлять дополнительными параметрами метода сеток.Однако пользоваться ею намного сложнее, чем функцией P d e s o l v e O , и поэтому в нашей книге мы не будем на ней особо останавливаться.Глава 13. Уравнения в частных производных339В качестве примера использования этой новой функции Mathcad 11 (листинг 13.4) используем то же самое одномерное уравнение теплопроводности(5) с граничными и начальными условиями (6) и (7).| Листинг 13.4.
Решение одномерного уравнения теплопроводности|D:=0.1L:=lT:=10Givenut(x, t) • D - u x x ( x , t)u(x,0)= Ф { x - 0.45 )u ( 0 , t) • 0-ф ( х - 0.55 )u ( L , t) = 0M:w:u := P d e s o l v e | u . x , [i.t.l1,100,10TДля корректного использования функции pdesolve предварительно, послеключевого слова Given, следует записать само уравнение и граничные условия при помощи логических операторов (для их ввода в Mathcad существуетспециальная панель). Обратите внимание, что уравнение должно содержатьимя неизвестной функции u(x,t> вместе с именами аргументов (а не так,как она записывается в пределах встроенной функции Pdesolve).
Для идентификации частных производных в пределах вычислительного блока следуетиспользовать нижние индексы, например, и ж ( , ь ) для обозначения второйпроизводной функции и по пространственной координате х..Как видно из рис. 13.14, на котором изображены результаты расчетов полистингу 13.4, встроенная функция с успехом справляется с уравнениемдиффузии, отыскивая уже хорошо знакомое нам решение.
Заметим, что использование встроенной функции pdesolve связано с довольно громоздкими вычислениями, которые могут отнимать существенное время.ПримечаниеКак Вы можете заметить, выбирать величину шага по пространственной и временной переменным может как сам алгоритм, так и пользователь (неявным образом, через число узлов сетки). Читателю предлагается повторить вычислениялистинга 13.4 для различных комбинаций параметров (главным образом, числаузлов сетки), чтобы проверить, в каких случаях алгоритм встроенной функциисправляется с задачей, выдавая верное решение, а в каких дает сбой.Приведем еще один пример применения функции Pdesolve для решенияуравнений в частных производных.
Рассмотрим одномерное волновое уравнение, которое описывает, например, свободные колебания струны музыкального инструмента:Эди(х, t)ъег'сЭ2и(х, t) .дх2(11)Часть III. Численные методы340Рис. 13.14. Решение уравнения диффузии теплапри помощи встроенной функции Pdesolve (листинг 13.4)Здесь неизвестная функция u(x,t) описывает динамику смещения профиляструны относительно невозмущенного (прямолинейного) положения, а параметр с характеризует материал, из которого изготовлена струна.Как Вы видите, уравнение (И) содержит производные второго порядка какпо пространственной координате, так и по времени. Для того чтобы можнобыло использовать встроенную функцию pdesolve, необходимо переписатьволновое уравнение в виде системы двух уравнений в частных производных,введя вторую неизвестную функцию v=ut.
Программа для решения волнового уравнения приведена в листинге 13.5, а результат — на рис. 13.15.Листинг 13.5. Решение волнового уравненияс:=1T:=lL := 2 • лGivenv t ( x , t) = с- u x x (х , t !u,-(x,t) в v ( x , t)u ( x , 0) s s i nufO, t) a 0v(x,0) =0u ( L , t) и 0Глава 13. Уравнения в частных производных:= Pdesolve0LN341°Т«(*;.Рис.
13.15. Решение волнового уравнения (листинг 13.5)13.3.2. Эллиптические уравненияРешение эллиптических уравнений в частных производных реализованотолько для единственного типа задач — двумерного уравнения Пуассона.Это уравнение содержит вторые производные функции и ( х , у ) по двум пространственным переменным:, у)э : и(х, у)— += -f(x, у)ду аЭх(12)Уравнение Пуассона описывает, например, распределение электростатического поля и(х,у) в двумерной области с плотностью заряда f ( x , y ) или(см. разд.
13.1.2) стационарное распределение температуры и(х,у) на плоскости, в которой имеются источники (или поглотители) тепла с интенсивностью f ( х , у ) .ПримечаниеjНесмотря на то, что применение встроенных функций, описанных в данномразделе, анонсировано разработчиками Mathcad только для уравнения Пуассона, их можно применять и для решения других уравнений, даже необязательноэллиптического типа. О том, как осуществить такие расчеты, написано вконце данного раздела.342Часть III. Численные методыУравнение Пуассона с нулевымиграничными условиямиКорректная постановка краевой задачи для уравнения Пуассона требует задания граничных условий.
В Mathcad решение ищется на плоской квадратной области, состоящей из (м+1)х(м+1) точек. Поэтому граничные условиядолжны быть определены пользователем для всех четырех сторон упомянутого квадрата. Самый простой вариант — это нулевые граничные условия,т. е. постоянная температура по всему периметру расчетной области. В таком случае можно использовать встроенную функцию muitigrid.П muitigrid(F,ncycie) — матрица решения уравнения Пуассона размера(м+1)х(м+и на квадратной области с нулевыми граничными условиями;• F — матрица размера (M+1)X(M+D, задающая правую часть уравненияПуассона;• ncycie — параметр численного алгоритма (количество циклов в пределах каждой итерации).Внимание!Сторона квадрата расчетной области должна включать точно М=2" шагов, т.е.2 п +1 узлов, где п — целое число.Параметр численного метода ncycie в большинстве случаев достаточновзять равным 2.
Листинг 13.6 содержит пример использования функцииmuitigrid для расчета краевой задачи на области ззхзз точки и точечнымисточником тепла в месте, задаваемом координатами (15,20) внутри этойобласти.Листинг 13.6. Решение уравнения Пуассонас нулевыми граничными условиямиМ := 32FM,M'•- 0Fl5,20 := Ю 4G — m u i t i g r i d ( - F , 2)В первой строке листинга задается значение м=32, в двух следующих строкахсоздается матрица правой части уравнения Пуассона, состоящая из всех нулевых элементов, за исключением одного, задающего расположение источника.
В последней строке матрице G присваивается результат действияфункции muitigrid. Обратите внимание, первый ее аргумент сопровождается знаком "минус", что соответствует записи правой части уравнения Пуас-Глава 13. Уравнения в частных производных343сона (11). Графики решения показаны на рис. 13.16 и 13.17 в виде трехмерной поверхности и линий уровня, соответственно.Уравнение Пуассона с произвольнымиграничными условиямиВ более сложных случаях, например для решения краевой задачи с ненулевыми условиями на границах, следует использовать другую встроеннуюфункцию relax, имеющуюся в Mathcad.П reiax(a,b,c,d,e,F,v,rjac} — матрица решения дифференциальногоуравнения в частных производных на квадратной области, полученногос помощью алгоритма релаксации для метода сеток;•a , b , c , d , e — квадратные матрицы коэффициентов разностной схемы,аппроксимирующей дифференциальное уравнение;• F — квадратная матрица, задающая правую часть дифференциальногоуравнения;• v — квадратная матрица граничных условий и начального приближения к решению;•rjac — параметр численного алгоритма (спектральный радиус итераций Якоби).Рис.