ПЗ - Численное модлирование распределения температурного поля

Введение
Изучение окружающего мира является одним из важнейших направлений для человечества. На смену менее точным методам, таким как наблюдение, математический расчет и экспериментальное подтверждение пришли современные методы математического моделирования, позволяющие проследить изменение состояния мира через время или смоделировать интересующий процесс, воссоздание которого в реальности очень сложно либо ограничено технически и финансово.
Несмотря на огромные возможности математического моделирования человек не способен в полной мере их использовать из-за естественных физических ограничений, таких как человеческая внимательность и концентрация, ограниченность ресурсов, например, бумага и чернила, а также ограниченность времени.
С задачей математического моделирования в современном мире хорошо справляются информационные технологии. С учетом стремительного роста мощностей компьютеров можно сказать, что скорость выполнения расчетов будет повышаться.
Основной проблемой для изучения мира путем математического моделирования с использованием информационных технологий является внутренняя компьютеров архитектура, которая не может производить сложные математические вычисления интегралов, дифференциальных уравнений, уравнений в частных производных, а также выполнять сложные математические преобразования функций.
////


1. ОБЗОР ЧИСЛЕННЫХ МЕТОДОВ МОДЕЛИРОВАНИЯ В МЕХАНИКЕ
1.1 Общая характеристика численных методов

Одним из доступных способов познания мира, в частности явлений, процессов и т.д., является математическое моделирование изучаемого объекта. Математическая модель является совокупностью уравнений, связывающих между собой параметры, характеризующие состояние исследуемого явления. В результате анализа составляется конкретная математическая задача, другими словами – набор уравнений, подлежащих решению.
За счет активного развития компьютерных технологий для решения анализа математических моделей прибегают к численным методам, позволяющим автоматизировать процесс.
Численные методы – это методы решения математических задач в численном виде, представление как исходных данных в задаче, так и её решения в виде числа или набора чисел.
Основами для вычислительных методов являются:
– решение систем линейных уравнений;
– интерполирование и приближённое вычисление функций;
– численное интегрирование;
– численное решение системы нелинейных уравнений;
– численное решение обыкновенных дифференциальных уравнений;
– численное решение уравнений в частных производных (уравнений математической физики);
– решение задач оптимизации.
Область математики, изучающая численные методы описывает основные алгоритмы построения моделей окружающего мира, а также способы его изучения, позволяя строить модели сложнейших физических, технических и биологических процессов.
Наблюдая за окружающим миром человек описывает их вербально, пытаясь понять суть явления, находя закономерности, после чего строиться математическая модель.
////


1.2 Погрешности в численных методах
Численные методы являются методами для поиска приближенного решения задач, а значит – результат вычислений всегда имеет некоторую погрешность и чем она меньше, тем точнее метод. Существует несколько видов погрешностей.
Погрешность модели. Современный мир слишком сложен для изучения во всем многообразии, во всей полноте его взаимосвязей, в том числе и незначительных. Любая наука, изучающая природу, не изучает ее напрямую, как было сказано ранее, изучает модели, созданные на основе процесса или явления. Простыми словами модель – это идеальное описание явления, не включающее в себя незначительные свойства и становится очевидно, что процесс моделирования сопровождается упрощением явления, а значит – добавлением погрешностей в результат описания модели.
Математическая модель создается на языке математики, но оценка погрешности математической модели есть прерогатива не математики, а той науки, в рамках которой изучается явление.
Погрешность исходных данных. Как правило, математическая модель содержит некоторые параметры, зависящие от исходных данных, которые определяются экспериментальным путем, а значит – являются не точными.
Погрешности в решении, обусловленные моделированием и исходными данными, называются неустранимыми. Они не зависят от математики и присутствуют, даже если решение поставленной математической задачи найдено точно.
////


1.3 Численные методы в механике
Механика – раздел физики, наука, изучающая движение материальных тел и взаимодействие между ними; при этом движением в механике называют изменение во времени взаимного положения тел или их частей в пространстве.
Численные методы в механике применяются во множестве областей:
1) аэромеханика;
2) вычислительная механика;
3) газовая и волновая динамика;
4) геометрия и топология;
5) гидромеханика;
6) механика композитов;
7) прикладная механика;
8) теоретическая механика и мехатроника;
9) теория пластичности;
10) теория упругости;
11) термодинамика;
12) и др.
/////


1.4 Методы решения численных методов в термодинамике
Для решения численных методов в термодинамике используется уравнение теплопроводности, которое, в общем виде, в пространстве с декартовыми координатами x = (x1, …, xn) имеет вид:

где a – коэффициент теплопроводности, f(x, t) – функция тепловых источников.
Уравнение теплопроводности называют однородным если f(x, t) = 0, т.е. система является теплоизолированной.
Для различных задач применяются различные модификации уравнения.

1.4.1 Задача Коши для уравнения теплопроводности. Одной из основных задач теории дифференциальных уравнений является задача Коши. Суть задачи состоит в нахождении интеграла уравнения, удовлетворяющего некоторым начальным условиям.
Задача Коши составляется при анализе процессов, изменяемых во времени, например – распределение тепла в плоскости. Этим обуславливается выбор терминологии: начальные значения и анализируемые значения, т.е. изменяющиеся во времени.
От краевых задач задача Коши отличается тем, что область, в которой определено искомое решение заранее не указывается. Тем не менее, задачу Коши можно рассматривать как одну из краевых задач.
////

Введение…………………………………………………………………. 3
1. ОБЗОР ЧИСЛЕННЫХ МЕТОДОВ МОДЕЛИРОВАНИЯ В МЕХАНИКЕ……………………………………………………………...
4
1.1 Общая характеристика численных методов…………………….. 4
1.2 Погрешности в численных методах……………………………… 6
1.3 Численные методы в механике……………………………………. 7
1.4 Методы решения численных методов в термодинамике………. 8
2. АЛГОРИТМИЧЕСКИЙ АНАЛИЗ ЗАДАЧИ………………………. 12
3. ПРОГРАММНАЯ РЕАЛИЗАЦИЯ ЗАДАЧИ……………………….. 17
3.1 Решение задачи с помощью пк MathCAD………………………… 17
3.2 Решение задачи с помощью пк ФОРТРАН………………………. 20
3.3 Решение задачи с помощью пк ANSYS…………………………… 24
Заключение………………………………………………………………. 31
Список использованных источников…………………………………. 32
Приложение А Код программы в среде MathCAD………………….. 33
Приложение Б Код программы в среде ФОРТРАН…………………. 35
////


2. АЛГОРИТМИЧЕСКИЙ АНАЛИЗ ЗАДАЧИ

Если требуется найти распределение температуры в твердом теле в случае, когда температура зависит от двух или трех пространственных координат, самый очевидный подход - попытаться получить точное решение основного уравнения. Уравнение стационарной теплопроводности в твердом теле с постоянным коэффициентом теплопроводности при отсутствии внутреннего тепловыделения имеет вид:

Это уравнение Лапласа. Уравнение Лапласа является линейным дифференциальным уравнением в частных производных. Известно несколько стандартных методов его решения. Один из них, например, метод разделения переменных, особенно полезен для решения задач теплообмена.
После того как тем или иным методом распределение температуры найдено, тепловой поток определяется с помощью закона Фурье. В двух- и трехмерных системах этот закон удобнее всего выразить в векторной форме:

где - градиент температуры (скалярной величины).
Градиент скалярной величины, такой, как температура, является вектором, т. е., согласно векторной записи закона Фурье, плотность теплового потока - это вектор. Обычно не рассматривается плотность теплового потока как вектор, поскольку она имеет размерность мощности на единицу площади, а ни одна из этих величии не является вектором. Однако удобно вообразить, что тепло «течет» в некотором направлении; поэтому величину часто называют вектором плотности теплового потока.
Вектор плотности теплового потока обладает важным геометрическим свойством, присущим градиентам: он направлен по нормали к изотерме, линии постоянной температуры, во всех точках твердого тела. Для иллюстрации этого свойства на рисунке 2 показаны несколько изотерм и типичных векторов плотности теплового потока в точках А, В и С твердого тела. Длина каждого из трех векторов плотности теплового потока пропорциональна местному градиенту температуры. Это значит, что в области тесного расположения изотерм градиент велик и плотность теплового потока также велика. В области, где расстояние между изотермами больше, плотность теплового потока соответственно меньше. На рисунке 1 плотность теплового потока в точке А больше, чем в точке В, где градиент температуры меньше.
/////


Решение: сначала нанесем на тело сетку с квадратными ячейками, как показано на рисунке 6. Пронумеруем узлы с неизвестными температурами от 1 до 16. Основная ячейка сетки представляет собой квадрат со стороной Δх = Δу = 250 мм. Температуры в узлах на 1-ой поверхности известны.


Рисунок 6 – Создание сетки

Система 16 разностных уравнений баланса энергии записывается следующим образом:
Узел 1:
(q∙∆x)/∆x+T_2-T_1=0;
Узел 2:
T_2=〖0,25T〗_1+〖0,25T〗_3+〖0,25T〗_7+250;
Узел 3:
T_3=〖0,25T〗_2+0,25T_4+〖0,25T〗_8+250;
Узел 4:
T_4=0,25T_3+〖0,25T〗_5+〖0,25T〗_9+250;
Узел 5:
(1000∙T_10)/2+T_4+(20∙h∙∆x)/k-T_5∙(2+(h∙∆x)/k)=0;
Узел 6:
(q∙∆x)/∆x+T_7-T_6=0;
Узел 7:
T_7=0,25T_2+〖0,25T〗_6+〖0,25T〗_8+〖0,25T〗_11;
Узел 8:
(T_9∙T_13)/2+T_3+T_7+(20∙h∙∆x)/k-T_8∙(3+(h∙∆x)/k)=0;
Узел 9:
(T_8∙T_10)/2+T_4+(20∙h∙∆x)/k-T_9∙(2+(h∙∆x)/k)=0;
Узел 10:
(T_5∙T_9)/2+T_4+(20∙h∙∆x)/k-T_10∙(2+(h∙∆x)/k)=0;
Узел 11:
(q∙∆x)/∆x+T_12-T_11=0;
Узел 12:
T_12=0,25T_7+〖0,25T〗_11+〖0,25T〗_13+〖0,25T〗_15;
Узел 13:
(T_8∙T_16)/2+T_12+(20∙h∙∆x)/k-T_13∙(2+(h∙∆x)/k)=0;
Узел 14:
T_14=0,0008T_11+0,0008T_15+27,955;
Узел 15:
(T_14∙T_16)/2+T_12+(20∙h∙∆x)/k-T_15∙(2+(h∙∆x)/k)=0;
Узел 16:
(T_13∙T_15)/2+(20∙h∙∆x)/k-T_16∙(1+(h∙∆x)/k)=0.

Полученную систему из 16 уравнений используем для определения неизвестный температур в программе MathCAD с помощью функции Find() (см. приложение А):
////


3.2 Решение задачи с помощью пк ФОРТРАН

Для численного решения задач на ЭВМ очень удобен итерационный метод, основанный на непосредственном определении температуры в каждом узле из разностного уравнения баланса энергии для этого узла. Например, если мы рассматриваем уравнение баланса энергии для внутреннего узла двумерного твердого тела, то получаем уравнение:

Разрешая это уравнение относительно температуры в узле 0, получаем:

Это соотношение типично для внутреннего узла в твердом теле с постоянными теплофизическими свойствами при отсутствии внутреннего тепловыделения, если применяется сетка с квадратными ячейками. Аналогичное соотношение получается для температуры в узле, расположенном на границе тела. Например, если узел 0 находится на границе, где происходит конвективный теплообмен с окружающей средой. температуру Т0 можно найти из уравнения:

Итак, температуру в каждом узле можно выразить через температуры в соседних узлах. Число полученных соотношений равно числу узлов с неизвестными температурами.
При использовании итерационного метода последовательно выполняются следующие четыре операции.
Операция 1. Выводит разностные уравнения, записав баланс энергии для каждого узла с неизвестной температурой. Из каждого уравнения выражают в явном виде температуру узла, для которого составлялся баланс энергии. Уравнения для всех внутренних узлов одинаковы по форме. Уравнения для граничных узлов будут различными в зависимости от типа граничных условий в конкретной задаче.
Операция 2. Задают ряд значений температур во всех узлах. Если задача будет решаться вручную, разумная начальная оценка всех температур позволит снизить затраты времени на вычисление истинных значений температуры в каждом узле. Если проводится численный расчет на ЭВМ, удобно принять все начальные температуры равными нулю.
Операция 3. Вычисляют новые значения температур, используя уравнения, полученные при операции 1. Как только получено новое значение какой-либо температуры, немедленно заменяют ее старое значение новым, так что новые значения температур в узлах все время вычисляют с использованием самого последнего приближения для остальных температур. Это позволяет уменьшить время сходимости решения к конечным стационарным значениям температур. Этот частный вид итерационного метода часто называют методом Гаусса—Зайдем.
/////


3.3 Решение задачи с помощью пк ANSYS

1. Создание плоскости
Preprocessor – Modeling – Create – Areas – Rectangle – By Dimensions


2. Создание материалов
Preprocessor – Material Props – Material Models/ Thermal –
Conductivity – Isotropic
Material – New Model… (Создание нового материала)

3. Создание конечного элемента
Preprocessor – Element Type – Add/ Add…/ Thermal Mass – Solid/
Quad 4 node 55

4. Разбитие плоскости на конечные элементы
Preprocessor – Meshing – MeshTool/ Set… (Выбор материала) / Mesh
(Создание сетки)


5. Приложение параметров
5.1 Задание конвекции
Thermal – Convection – On Lines
(К боковой поверхности стержня)
Вторая строка – коэффициент теплообмена
Четвертая строка – Температура внешней среды
////

Приложение Б
Код программы в среде Фортран

!Программа iter
program iter
Implicit None
real, dimension (16) :: T
real, dimension (16) :: TT
real summ, Dx, h, k
real TOLER
integer I
!Допустимое отклонение
TOLER = 1
!Назначаем на первом шаге расчёта все температуры равные 0
do I = 1, 16
T(I) = 0
enddo
!Пока не достигли допустимого отклонения делаем цикл

!summ = 20
Dx = 0.25
h = 500
k = 0.2


do while (summ > 16)

!Узел 1:
TT(1) = 800*Dx/k+T(2)
!Узел 2:
TT(2) = 0.25*T(1)+0.25*T(3)+0.25*T(7)+250
!Узел 3:
TT(3) = 0.25*T(2)+0.25*T(4)+0.25*T(8)+250
////

Задача: требуется найти стационарное распределение температуры от всех поверхностей двумерного тела, показанного на рисунке 5.

Решение: сначала нанесем на тело сетку с квадратными ячейками, как показано на рисунке 6. Пронумеруем узлы с неизвестными температурами от 1 до 16. Основная ячейка сетки представляет собой квадрат со стороной Δх = Δу = 250 мм. Температуры в узлах на 1-ой поверхности известны.

1. Токочаков, В.И. Моделирование, оптимизация и управление теплотехническими системами: лабораторный практикум по одноименному курсу для студентов специальности 1-43 01 05 ≪Промышленная теплоэнергетика≫ дневной и заочной форм обучения: в 3 ч. Ч. 2 / В.И. Токочаков. – Гомель: ГГТУ им. П. О. Сухого, 2009. – 41 с.
2. Цветков Ф.Ф. Тепломассообмен: учебное пособие для вузов / Ф.Ф. Цветков, Б.А. Григорьев – Москва: Издательский дом МЭИ, 2005. – 550 с.
3. Бахвалов Н.С. Численные методы. Учебник / Н.С. Бахвалов, Н.П. Жидков, Г.М. Кобельков – БИНОМ. Лаборатория знаний, 2012. – 636 с.
4. Зенков А.В. Численные методы: учебное пособие / А.В.Зенков. – Екатеринбург: Изд-во Урал. ун-та, 2016. – 124 с.
5. Рихтер Д. CLR via C#. Программирование на платформе Microsoft .NET Framework 4.5 на языке C# /Д. Рихтер. – С.Питербург 2012. – 896 c.
6. Уравнения математической физики [Электронный ресурс].: Научная библиотека. – Электронные данные. – Режим доступа.: http://scask.ru/q_book_emp.php?id=37/ – Дата доступа: 04.02.2020.