Содержание к диссертации
Введение
Глава 1. Методы и алгоритмы исследования оптимизационной модели распределения плотности источников тепла в однородной неподвижной среде 22
1.1. Построение и преобразование оптимизационной модели. Конечномерная аппроксимация и её регулярность по функционалу 22
1.2. Оптимизационная модель в одномерном случае. Алгоритм нахождения обменной матрицы с помощью функции Грина и численными методами 29
1.3. Двумерная модель. Алгоритм построения конечно-разностной схемы для решения прямой задачи теплопроводности 35
1.4. Алгоритм построения конечно-разностной схемы для решения прямой задачи теплопроводности в трёхмерном случае 39
1.5. Численное исследование оптимизационной модели. Результаты вычислительных экспериментов 41
Глава 2. Методы и алгоритмы исследования оптимизационной модели распределения плотности источников тепла в неодно родной неподвижной среде 50
2.1. Оптимизационная модель, её конечномерная аппроксимация и регулярность по функционалу 50
2.2. Одномерная модель. Алгоритм численного решения прямой задачи теплопроводности на отрезке 52
2.3. Двумерная модель 59
2.4. Трёхмерная модель 61
2.5. Результаты вычислительных экспериментов 62
Глава 3. Методы и алгоритмы исследования оптимизационной модели распределения плотности источников тепла в одно родной движущейся среде 68
3.1. Формулировка и преобразование краевой задачи установившегося теплообмена 69
3.2. Построение оптимизационной модели распределения плотности источников тепла в движущейся среде 73
3.3. Конечномерная аппроксимация оптимизационной модели 76
3.4. Регулярность конечномерной аппроксимации по функционалу 77
3.5. Алгоритм построения конечно-разностной схемы для вычисления поля скоростей 82
3.6. Алгоритм нахождения элементов обменной матрицы aij 85
3.7. Численное исследование оптимизационной модели. Результаты вычислительных экспериментов 86
Глава 4. Описание программно-инструментального комплекса HeatCore 93
4.1. Технические характеристики HeatCore 93
4.2. Используемые обозначения 96
4.3. Основы работы в HeatCore 97
Заключение 111
Литература
- Двумерная модель. Алгоритм построения конечно-разностной схемы для решения прямой задачи теплопроводности
- Одномерная модель. Алгоритм численного решения прямой задачи теплопроводности на отрезке
- Результаты вычислительных экспериментов
- Алгоритм построения конечно-разностной схемы для вычисления поля скоростей
Двумерная модель. Алгоритм построения конечно-разностной схемы для решения прямой задачи теплопроводности
Работа носит теоретический характер. В ней поставлена новая задача нахождения оптимального распределения плотности источников тепла и разработан метод её решения. В перспективе, предложенные алгоритмы можно использовать при разработке расширений для CAE-систем. На практике программный комплекс может быть использован для синтеза оптимальных систем обогрева жилых и производственных помещений (теплиц и т.д.).
Область исследования. Содержание диссертации соответствует паспорту специальности 05.13.18 «Математическое моделирование, численные методы и комплексы программ» по следующим областям исследований:
Развитие качественных и приближённых аналитических методов исследования математических моделей.
Разработка, обоснование и тестирование эффективных вычислительных методов с применением современных компьютерных технологий.
Реализация эффективных численных методов и алгоритмов в виде комплексов проблемно-ориентированных программ для проведения вычислительного эксперимента.
Апробация работы. Результаты работы докладывались на следующих научных конференциях: XXIV Международная научная конференция — «Математические методы в технике и технологиях (ММТТ-24)» (Национальный технический ун-т Украины «КПИ», г. Киев, 2011 г.); Международная конференция «Комплексный анализ и его приложения в дифференциальных уравнениях и теории чисел» (ИПК НИУ «БелГУ», г. Белгород, 2011 г.); Международная молодёжная конференция «Прикладная математика, управление и информатика» (ИПК НИУ «БелГУ», г. Белгород, 2012 г.); XXV Международная научная конференция — «Математические методы в технике и технологиях (ММТТ-25)» (Национальный технический ун-т «ХПИ», г. Харьков, 2012 г.); Международная конференция «Дифференциальные уравнения и их приложения» (ИПК НИУ «БелГУ», г. Белгород, 2013 г.); VIII Международная научно-практическая конференция «Инновационное развитие: физико-математические и технические науки» (г. Москва, август 2014); а также на семинарах Орловского государственного университета, Белгородского государственного технологического университета им. В.Г. Шухова, НИУ «БелГУ».
Публикации. Основные результаты диссертации отражены в 9 публикациях, 3 из которых в изданиях, рекомендованных ВАК. В результате работы над диссертацией создан ряд программ для ЭВМ, три из которых получили государственную регистрацию.
Работа использует модель стационарного теплообмена [71, 109], основанную на ряде предположений, определяющих условия, при которых решается наша задача. Существует три механизма распространения тепла: теплопередача посредством межмолекулярных взаимодействий (теплопроводность), конвекция и излучение. Труднее всего моделируется свободная конвекция, возникающая в гравитационном поле. В работе рассматриваются только стационарные тепловые процессы, в которых свободной конвекцией, как правило, можно пренебречь. Теплообмен электромагнитным излучением играет существенную роль при температурах, превышающих 500 C, и мы его также не учитываем. Теплопроводность описывается законом Фурье q= —кЧТ, связывающим вектор плотности потока тепла q с градиентом температуры Т. Здесь к — коэффициент теплопроводности, являющийся характеристикой вещества, заполняющего некоторую область. Из закона Фурье следует (см., например, [57, с. 7]), что температурное поле Т(х) в области, заполненной неподвижной средой, удовлетворяет стационарному уравнению теплопроводности V {кЧТ) — f(x) = 0, где f(x) — плотность мощности источников тепла, находящихся в области. Аналогичное уравнение справедливо при наличии вынужденной конвекции среды (см. ниже параграф 3.1).
Теплообмен с окружающей средой описывается законом Ньютона-Рихмана Q = р(х)(Т(х) — То)AS [45, с. 26], который выражает величину теплового потока Q через элемент границы (стенки) области площадью AS около точки х с разностью между температурой стенки Т(х) и температурой окружающей среды То. Величину /І(ж) называют коэффициентом теплоотдачи в окружающую среду. Закон Ньютона-Рихмана совместно с законом Фурье и предположением, что температура стенки совпадает с температурой внутренней среды у этой стенки, приводят к следующему краевому усло ЭТ т/ вию % Ь OL[X){T — То) = 0 для температурного поля Т(х), справедливо му на границе области. Здесь х коэффициент температуропроводности, X = х/(рс), р — плотность среды, с — её удельная теплоемкость. Величину а(х) мы в дальнейшем называем коэффициентом теплопередачи через границу [103]. В случае, когда среда внутри области жидкая или газообразная, тепловое взаимодействие с твердой стенкой является довольно сложным. В частности, предположение о том, что температура стенки совпадает с температурой внутренней среды у этой стенки нельзя считать справедливым. Тем не менее, краевое условие можно считать выполненным при соответствующем выборе коэффициента теплопередачи а(х). Все трудности в описании теплоотдачи в окружающую среду можно перенести на определение коэффициента теплопередачи а(х) (см., например, [96, с. 40]). В дальнейшем мы считаем, что температурное поле Т(х) внутренней среды удовлетворяет стационарному уравнению теплопроводности и краевому условию третьего рода, приведенному выше.
Одномерная модель. Алгоритм численного решения прямой задачи теплопроводности на отрезке
Таким образом, для любого сколь угодно малого а 0 найдется такой номер Na, что при п А , выполнено неравенство (Jn)min (7 + ск)(1 — ЦСЦС/ о" а)-1, которое и доказывает теорему.
Эта теорема точно также доказывается для других, упомянутых выше, односторонних модификаций рассматриваемой задачи.
Замечание 2. Регулярность по функционалу последовательности двусторонних задач Zo(n) легко устанавливается по приведенной схеме, если для задач Zo(n) и Z(n) удаётся доказать справедливость леммы 2. Её доказательство сводится к установлению существования точки минимума задачи Zo(n), для которой правые неравенства в (5) являются строгими. В пользу существования такой точки минимума говорит замечание 1. 1.2. Оптимизационная модель в одномерном случае. Алгоритм нахождения обменной матрицы с помощью функции Грина и численными методами коэффициенты теплопередачи соответственно через левую и правую границы, м/с. Решение данной краевой задачи может быть получено в аналитическом виде с использованием функции Грина [91]. Данная функция G(x,) для нашей задачи подчиняется следующим условиям (см. [77, с. 288]):
Воспользуемся этими свойствами для построения функции 6г(ж,) в явном виде. Графически, функция Грина и её производная имеют вид кусочно-линейных функций (рис. 1.1), если зафиксировать второй аргумент константой = о. Будем искать решение в виде: где Мі (ж), г 2(ж) являются однородными уравнениями, удовлетворяющие условиям и [ = 0, v,2 = 0. В общем виде их решения можно представить как щ{х) = А\х + В\, іі2(х) = А2Х + B i. Подставив функцию G(x,) в условие (1.9), получим:
Аналогично найдём функцию v,2(x) для правой границы отрезка. Очевидно, она должна удовлетворять условию (и2(х) + / (ж)) _,= 0. Опустив подробный вывод, запишем сразу + $ф U
А\(х + h)zi(), а х ; (ж — /2) 2( )5 « ж 6. Теперь найдём zi и Z2 такие, чтобы удовлетворить требованиям (1.8), (1.10). В точке х = функция Грина непрерывна, а производная терпит разрыв, равный 1, значит
Выразим из второго уравнения A2Z2 и подставим в первое: Приступим к определению элементов обменной матрицы CLij = (Gej,ei). Gej можно вычислять аналитически, а не интегрировать численно, тем самым существенно уменьшить в дальнейшем процессорное время при расчёте:
Изложим также численный метод построения обменной матрицы а . Это необходимо для дальнейшего сравнения алгоритмов в ходе получения элементов aij и проверки применимости численного способа, так как в многомерном случае решение краевой задачи (1.2) возможно только с использованием численных методов. Решение (1.6), (1.7) будем находить с использованием конечно-разностного метода [74, с. 93].
Аппроксимируем функцию и сеточной функцией yk (к = 0,..., N) внутри (О к N) отрезка [a, Ь] трёхточечным оператором: около узла — J f(x)dx, тогда разностная схема обладает свойством консер Хк-1/2 ъ вативности и гарантирует, что а\уо + OLZVN = f f{x)dx. Из полученных раз а ностных уравнений формируется система линейных уравнений с неизвестны ми ук, которую можно решать с использованием итерационных методов или метода Гаусса. На рисунке (1.2) представлена блок-схема численного решения задачи (1.6), (1.7). Решение данной задачи с / = е , по сути, и представляет собой у = Gej. Какой бы способ из описанных выше двух не использовался для получения Gej, для окончательного построения обменной матрицы а нужно произвести ещё одно интегрирование а = J Gejdx.
Подобласти также выберем разбиением отрезка [а, Ь] на одинаковые части таким образом, чтобы N было кратным п (N = 2п, Зп,4п,...). На самом деле, сначала выбирается п, а затем N такое, чтобы на границах Dj обязательно были узлы разностной схемы.
Блок-схема решения прямой задачи теплопроводности на отрезке 1.3. Двумерная модель. Алгоритм построения конечно-разностной схемы для решения прямой задачи теплопроводности
Поскольку построение функции Грина в многомерном случае затруднительно, решение (1.2) будем находить численно с помощью метода конечных разностей [34, 42]. задают число разбиений области D. В наших экспериментах Пі,П2 будут сравнительно небольшими числами, меньшими 100, ибо при сильном разбиении очень большим является в дальнейшем при решении обратной задачи время работы симплекс-метода. На рис. 1.3 показан пример такого разбиения.
Результаты вычислительных экспериментов
Введём функцию и(х) = Т(х) — То, где Т(х) — установившаяся температура в точке х ограниченной связной области D, а То — температура окружающей среды. Эта функция должна удовлетворять уравнению где х коэффициент температуропроводности среды, v(x) — поле скоростей среды, которое предполагается подчиненным условию divz? = О, а f(x) — плотность источников тепла в области D. Для формулировки краевых условий разобьём границу области D на три части 3D = Г+иГоиГ_, где Г+ — часть границы, которая является входом среды в область D, Г_ — часть границы, являющейся выходом (стоком) среды, а Го — часть непроницаемой границы. Справедливы следующие соотношения: единичный вектор внешней нормали к границе 3D. Последнее из этих трех условий означает, что в область D поступает вещество внешней среды с температурой То. В дальнейшем мы предполагаем, что в нашем процессе присутствует поток тепла через непроницаемую для среды границу, равный а(х) и(х) (х Є Го), где а(х) 0 коэффициент теплопередачи через Го. Следующие краевые условия выражают тепловой баланс области D с окружающей средой.
Если полученная краевая задача имеет единственное решение, то её можно считать удовлетворительной моделью процесса установившегося теплообмена при наличии конвекции.
Однако вопрос о существовании и единственности решения краевой задачи (3.1), (3.2) сложен для анализа. В частном случае, когда поле скоростей среды является потенциальным, эту задачу можно преобразовать к более простому виду.
Оператор Lw = — xAif + ( \7(p2/(4%))ii;, действующий в пространстве L,2(D) на достаточно гладкие функции w{x), подчиненные краевым условиям (3.5), является самосопряжённым и положительно определённым, а поэтому имеет ограниченный обратный оператор, определённый на L2(D). Отсюда следует, что краевая задача (3.3), (3.5) имеет единственное решение, которое можно выразить через функцию Грина оператора L.
В оптимизационной задаче (3.12) S(D) является классом функций, среди которых разыскивается решение д(х). Этот класс желательно выбирать как можно шире, поскольку при этом можно рассчитывать на уменьшение экстремального значения J{g}. Самым широким из возможных классов S(D) можно было бы считать пространство суммируемых функций Li(D), поскольку на нем естественно определяется функционал J{g}. Однако на этом пространстве не определён оператор G и задача (3.12) становится сформулированной не корректно. Корректной эта формулировка будет, если в качестве S(D) взять пространство функций интегрируемых с квадратом L2(D). Ниже мы считаем, что S(D) С L,2(D). В особо оговоренных случаях считается, что S(D) = L2(D).
Если найдено точное или приближенное решение до(х) задачи (3.12), то разыскиваемая оптимальная плотность источников тепла равна /min(ж) = е х х до(х). (3.13) Замечание. Возможны и другие формулировки оптимизационной задачи. Иногда ограничения сверху в (3.12) отсутствуют или не существенны. Такую оптимизационную задачу назовем односторонней. В задаче (3.12) возможно появление дополнительных условий, которые могут привести к некоторым новым её модификациям. Одна из естественных модификаций состоит в требовании невозможности расположения источников тепла в некоторой части области D, т.е. возникает дополнительное требование: д(х) = 0 при х Є Do С D. Такую модификацию назовем задачей с ограничениями на локализацию источников. Еще одна модификация связана с присутствием некоторых фиксированных источников до оптимизации. При этом функция д{х) в (3.12) состоит из двух слагаемых, одно из которых известная функция, а второе подлежит определению. Такую задачу назовем задачей с фиксированными источниками. По форме она мало отличается от задачи (3.12). К функционалу J{g} добавляется постоянное слагаемое, а функции #і(ж), #2(ж) изменяются на однозначно определяемые слагаемые. Возможны различные комбинации рассмотренных выше модификаций задачи (3.12).
Мы не обсуждаем здесь вопрос о существовании точного решения задачи (3.12) и её модификаций. Вообще говоря, подобные задачи не обладают свойством единственности решения (см. гл. 1, замечание 1) Ответы на эти вопросы сильно зависят от выбора класса функций S(D). Для практических нужд во многих случаях достаточно научиться находить так называемое квазирешение
Алгоритм построения конечно-разностной схемы для вычисления поля скоростей
Задание данных краевых условий осуществляется с использованием логических конструкций іґ( условие ). Краевое условие для нахождения функции температуры можно задать, например, следующим образом (рис. 4.5).
Данный пример снабжён комментариями в стиле C++. Краевое условие для нахождения поля скоростей задаётся в таком же диалоговом окне, но при этом необходимо следить, чтобы количество втекающего и вытекающего вещества в единицу времени было одинаковым. Если нужно задавать выход вещества, то переменной value необходимо присваивать отрицательное значение, а если вход, то положительное. Задание краевого условия в 2-мерном
Система следит за выполнением условия баланса протекающего вещества. Если баланс не выполняется, выдаётся предупреждение (рис. 4.7).
При создании нового файла расчёта всем входным параметрам присваиваются некоторые первоначальные значения. Количество разбиений N можно брать от десятков до нескольких тысяч. При задании профилей температур необходимо следить, чтобы значения функции минимального профиля m были меньше значений максимального профиля М, иначе алгоритм не найдет решение. Расчёт моделей для неоднородной среды отличается только тем, что в таблице параметров температуропроводность х задаётся в виде функции, а не как коэффициент. После задания всех входных данных расчёт производится нажатием клавиши F5 или выбором пункта главного меню Расчёт — Старт. При необходимости можно завершить расчёт нажатием F4 или выбором Расчёт — Остановить. Результат расчёта представляется в виде функций на декартовой плоскости (рис. 4.7). Красная и зеленая кривые - функции Мит соответственно. Оптимальная плотность источников имеет вид столбиков чёрного цвета, а результирующая температура — жёлтый цвет. Значение суммарной плотности источников тепла и интенсивность каждого источника показываются после расчёта в окне Вывод.
При расчёте двумерной области рекомендуется задавать количество источников не более 2000. Краевое условие а для решения прямой задачи теплопередачи задаётся в виде функции, определённой на границе. Расчёт оптимальной плотности источников в неоднородной среде отличается заданием переменной функции температуропроводности х. Для движущейся среды также необходимо задавать краевое условие для нахождения поля скоростей. Разбиения N{v}x, N{v}y можно задавать большими значениями — до 100.
Оптимальная плотность источников после расчёта представляется в виде столбиков (рис. 4.8), высота которых зависит от интенсивности источников. Источники тепла имеют оранжевый цвет, стоки - синий. Результирующая
Окно расчёта одномерной модели температура показывается в виде поверхности. В окно «Вывод» выводятся численные значения интегралов, основным из которых является первый, который обозначает суммарную интенсивность всех источников.
Чтобы увидеть поле скоростей (рис. 4.9) необходимо выбрать пункт меню Действия — Показать тепловую карту или нажать клавишу F10. Трёхмерная среда Рекомендованные числа разбиений на источники Nx, Ny, Nz — не более 11. Для расчёта поля скоростей можно задавать числа разбиений N{v}x, N{v}y, N{v}z до 30. Краевые условия и функции задаются в редакторе математических выражений.
В данном окне (рис. 4.9) синими линиями и летящими красными частицами показаны линии движения вещества. Менее нагретые области показаны зелёным цветом, более нагретые — красным. Всплывающая подсказка содер 107
Расчёт 2-мерной модели с движущейся средой жит информацию о параметрах среды в точке, в которой находится указатель курсора мыши. В неё выводятся текущие координаты, температура и плотность распределения источников.
В графическое окно (рис. 4.10) выводится информация о температуре внутри области и на границе, а также плотность распределения источников. Поскольку вся информация выводится в одно окно предусмотрено управление объектами сцены.
