|
Вернуться на страницу <Методические разработки>
Архив разработки  (110 Кб, MWord)
Решение задач теплопроводности и диффузии в пакете Mathematica.
Рассмотрим пример решения одномерного уравнения теплопроводности (диффузии) в системе Mathematica3.0. Как известно, такое уравнение в общем виде можно записать как
 .
Чтобы рассматриваемая задача была определена полностью, обычно добавляют еще начальные и граничные условия.
Пусть соответствующая задача формулируется следующим образом. Дан тонкий однородный стержень длиной l с теплоизолированной боковой поверхностью, начальная температура которого и(х,0)= Ах/l для 0<x<l. На конце стержня х = 0 температура поддерживается равной 0, а температура конца стержня x = l изменяется по закону u(l, t) = Ae -t, A=const.
Для простоты примем А = 1, l = 1. Тогда математическая формулировка данной задачи имеет вид
 Рис. 1
На рис.1 мы видим решение данной задачи в пакете Mathematica3.0. Названия стандартных функций вводятся с большой буквы, параметры функции – в квадратных скобках.
Построим график точного решения (формальные параметры вводятся с нижним подчеркиванием), которое имеет вид
 .
Соответствующий график с фрагментом вычислений показан на рис.2. Сумма ряда ограничена 25-ю членами. Очевидно, что численное и аналитическое решение идентичны.
 Рис. 2
Из рис.1 и 2 видны основные особенности полученного решения, в частности, быстрое остывание стержня до нулевой температуры и переход системы в стационарный режим.
Недостатком решения задач диффузии и теплопроводности в системе Mathematica является ограничение одномерными задачами. Двумерные задачи, даже стационарные, данная система не решает. Кроме того, система не решает даже относительно простые задачи, если в их условии или решении имеется сингулярность. Численное решение таких задач под силу системе MATLAB.
Визуализация решения задач математической физики в пакете MATLAB.
Для построения трехмерных графиков решения ДУ в частных производных в вычислительной системе MatLab имеется очень удобная система pdetool. Вызвать окно этой системы можно, набрав данную опцию или задав в командном окне область решения конкретной задачи командой
pderect([x1 x2 y1 y2]) - построение прямоугольника {x1 xx2, y1yy2},
pdecirc(x0,y0,R) - построение круга с центром (x0,y0) и радиусом R.
После набора в командном окне одной из данных опций на экране появляется окно, меню которого показано на рис.3.
Рис. 3
В качестве примера работы данной системы рассмотрим решение следующей задачи. Найти решение уравнения Лапласа Du = uхх + uуу = 0 в прямоугольнике D = {(x,y)| 0xa, 0yb}, удовлетворяющее следующим краевым условиям:
 .
Для определенности выберем значение параметров а = b = a = 1.
Опцией pderect задаем прямоугольник с координатами ([0 1 0 1]). Открывается окно PDE Toolbox, в котором заданный треугольник выделен серым цветом. Можно также задавать область решения непосредственно в данном окне в виде прямоугольника, круга, эллипса и многоугольника, выбирая соответствующие кнопки слева в нижнем ряду панели меню. В нижней части окна при различных действиях появляется краткая информация о данном действии.
MatLab решает двумерные задачи, используя метод конечных элементов. При этом расчетная область (в данном случае прямоугольник) разбивается на конечное число элементарных подмножеств стандартной формы (которые и называются конечными элементами). В системе MatLab конечные элементы представляют собой криволинейные треугольники. Чтобы задать первоначальное разбиение области решения на конечные элементы в окне PDE Toolbox достаточно в меню окна нажать кнопку D . Сделать исходное разбиение более мелким можно при помощи кнопки, стоящей рядом справа. После задания исходного разбиения переходим к заданию граничных условий с помощью кнопки ¶W. При этом границы области решения выделяются однонаправленными стрелками красного цвета. Щелчком левой кнопки мыши на каждой из стрелок открываем окно задания граничных условий Boundary Condition. Если на данном отрезке границы заданы условия Дирихле, этот отрезок выделяется красным цветом, условия Ньюмена - синим цветом, смешанные условия - зеленым цветом. Таким образом, после задания граничных условий мы по цвету граничных отрезков сразу можем определить тип условия, заданного на конкретном отрезке границы. Окно задания граничных условий показано на рис.4. Для нашей задачи на левой и правой границах прямоугольника заданы условия Дирихле с h = 1 и r = 1 и у, соответственно. На нижней и верхней границах заданы условия Ньюмана с g = q = 0. В верхней части окна при этом отображается вид граничного условия.
Рис. 4
После задания граничных условий задаем вид ДУ с помощью кнопки PDE. В открывшемся окне PDESpecification в данном случае выбираем опцию Еlliptic, а значения коэффициентов - с = 1, a = f = 0. Закрыв окно, нажимаем кнопку решателя =. Полученное решение отображается в исходном прямоугольнике в виде различно окрашенных областей. Цвет и интенсивность окраски соответствуют различным числовым значениям полученной функции.
Чтобы увидеть трехмерный график полученного решения, нажимаем кнопку справа от кнопки решателя. Открывается новое окно с трехмерным графиком полученного решения, который показан на рис.5.
 Рис. 5
Зафиксировав стрелку курсора на области графика и нажав левую кнопку мыши, полученную трехмерную фигуру можно вращать, добиваясь положения, при котором свойства решения видны наиболее отчетливо.
Граничные условия при использовании данной опции могут быть представлены в виде функции переменных х и у, однако если эта функция задана в виде степенного ряда с нелинейными членами, задание граничных условий становится более сложным. Это связано с тем, что все данные MatLab воспринимает как матрицы, поэтому возведение матрицы-столбца в степень воспринимается как некорректная задача. Чтобы этого избежать, необходимо степенные члены выражения задавать с точкой, что означает поэлементное действие.
В качестве примера рассмотрим предыдущую задачу с граничными условиями вида
 .
При решении нелинейное граничное условие задаем в виде y-y.^2. График полученного решения показан на рис.6.
С помощью пакета MatLab можно решать сингулярные задачи. Пример такой задачи можно увидеть, набрав в командном окне pdedemo7. В данном примере решается уравнение Пуассона  .
Точное решение данного уравнения  .
 Рис. 6
Программа определяет область решения в виде круга единичного радиуса с нулевыми условиями на границе и точечным источником в центре круга. Точечный источник задан в виде d -функции, которая аппроксимируется следующим образом: в пределах криволинейного треугольника, содержащего начало координат, функция принимается равной 1/S, где S - площадь треугольника; во всех остальных фрагментах функция принимается равной нулю.
Практическое задание.
1. Решить численно следующую задачу теплопроводности в системе Mathеmatica, построить трехмерный график полученного решения и сравнить его с графиком точного решения (аппроксимировать бесконечный ряд его конечной суммой).
На концах однородного изотропного стержня длиной l поддерживается нулевая температура. Предполагая, что стенки стержня теплоизолированы от окружающей среды, найти закон распределения температуры в стержне, если известно, что в начальный момент имелось распределение температуры
 .
Коэффициент диффузии равен D0.
l = n, u0 = 0.2n, D0 = 0.5n, tmax =1.5, где п – номер варианта.
Точное решение:
 .
2. Получить график решения задачи колебаний в системе MatLab.
Закрепленной по контуру однородной квадратной мембране со стороной l придали форму u(x, y,0) = sin(p x/l)sin(2p y/l). Начальная скорость точек мембраны постоянна и равна a/l, где а – постоянная, входящая в уравнение колебаний.
l = n, а = 0.5n, t = n, где п – номер варианта.
Обратить внимание на правильный ввод начального условия (с точкой после первого синуса и правильного значения параметра колебаний) для формы мембраны (опция Parameters в меню Solve). Вывести трехмерный график и наблюдать анимацию, установив соответствующий флажок в окне PlotSelection.
В системе Mathеmatica построить график точного решения (аппроксимировать бесконечный ряд его конечной суммой) для различных значений t (используя, например, опцию Do) и наблюдать анимацию полученного решения, выделив общую для всех графиков скобку и нажав Ctrl-Y. Точное решение:
3. Найти стационарное распределение температуры (построить график в системе MatLab) в прямоугольнике D={(x,y)| 0xa, 0xb}, если на границе прямоугольника поддерживается заданная температура:
u(x,0) = u(x,b) = 0, xI (0,a),
u(0, y) = y(b-y), u(a, y) = sin(p y/b), yI (0,b).
В системе Mathеmatica построить график точного решения (аппроксимировать бесконечный ряд его конечной суммой):
a = n, b = n/2, где п – номер варианта.
4. В MatLab рассмотреть пример решения сингулярной задачи эллиптического типа в pdedemo7.
|
