Содержание к диссертации
Введение
1. Анализ проблемы математического моделирования свободной конвекции
1.1. Математическое описание свободной конвекции 9
1.2. Формализованная запись основных уравнений 13
1.3. Типы постановок граничных условий на смоченной и свободной поверхностях
1.4. Экспериментальный подход определения гидротермических 19 характеристик
1.5. Аналитические точные и приближенные методы. 23
1.6. Численное интегрирование 30
1.7. Выводы, цель и задачи исследования 39
2. Приближенное аналитическое решение нестационарной задачи стокса .
2.1. Формулировка задачи и вывод основных уравнений 41
2.2. Решение задачи методом конечного интегрального синус преобразования Фурье
2.3. Анализ решения 52
2.4. Решение задачи в прямоугольнике 57
2.5. Стационарная постановка 59
2.6. Выводы 61
3. Численное интегрирование неоднородного бигармонического уравнения в стационарной и нестационарной постановках
3.1. Стационарная постановка 62
3.1.1. Постановка задачи и построение численной схемы 62
3.1.2. Вычисление нормы оператора перехода оптимального итерационного шага
3.1.3. Вычисление погрешности аппроксимации уравнения разностной схемой
3.1.4. Устойчивость и сходимость 72
3.2. Нестационарная постановка. 75
3.2.1. Постановка задачи и построение численной схемы 75
3.2.2. Вычисление нормы оператора перехода и 78 оптимального временного шага
3.2.3. Вычисление погрешности аппроксимации уравнения 82 разностной схемой
3.2.4. Устойчивость и сходимость 83
3.3. Улучшение решения на границе 84
3.4. Решение задачи в прямоугольнике 87
3.5. Выводы 92 4. Сравнение аналитического и численного подходов. 94
4.1. Стационарная постановка 94
4.2. Нестационарная постановка 96
4.3. Выводы 97
Заключение 99
Литература 100
Приложение 110
- Типы постановок граничных условий на смоченной и свободной поверхностях
- Численное интегрирование
- Решение задачи в прямоугольнике
- Вычисление нормы оператора перехода оптимального итерационного шага
Типы постановок граничных условий на смоченной и свободной поверхностях
Гравитационной конвекцией называется гидродинамическое явление, возникающее в поле всемирного тяготения, в связи с неодинаковой плотностью разных частей жидкости. Гидродинамическая сторона процесса проявляется в том, что частицы жидкости движутся не в пустом пространстве, а среди других подобных частиц, поэтому каждая частица при своем движении занимает место каких-либо других частиц, отодвинутых ею в сторону [1].
Причиной различия в плотности обычно являются различия в температуре или в составе, в частности в концентрации растворенных в жидкости примесей. Наиболее изученной формой гравитационной конвекции является температурная (тепловая) конвекция. Конвекция называется свободной, если те напряжения (в том числе нормальное давление), которые испытывает жидкость на ее границах, не совершают механической работы, т.е. если все границы жидкости неподвижны.
Первым исследователем, который научно подошел к явлениям тепловой конвекции был М. В. Ломоносов. Он первым правильно объяснил основной механизм метеорологических явлений. Систематическое исследование свободноконвективных течений началось с работ Дж. Томсона (1888), Б. Бенара (1900) и Л. Рэлея (1916), а математическая формализация явлений переноса, сопровождающих такие течения, впервые описана независимо друг от друга Обербеком (1879) и Буссинеском (1903) [2], которые установили, что основой такого описания служат фундаментальные уравнения Навье-Стокса в напряжениях, конкретизируемые реологическими уравнениями состояния, связывающими напряжения и тензор скоростей деформации. Наибольший практический интерес представляют так называемые ньютоновские жидкости, реологическое уравнение которых линейным образом связывает напряжение и скорость деформации среды, причем коэффициент пропорциональности /л называют динамической вязкостью. В этом случае уравнения Навье-Стокса без учета объемной вязкости в силу, как правило, небольших скоростей при свободноконвективном движении принимают в векторном виде классическую запись [3]: дт х где v, g - векторы скорости и ускорения силы тяжести; р - плотность среды; р - давление; т - текущее время; V, А - дифференциальные операторы "набла" и Лапласа. Уравнения (1.1) и (1.2) дополняются уравнением переноса теплоты, но без учета диссипации энергии опять же по причине небольших скоростей свободноконвективного течения
Во многих случаях представляет интерес исследование конвекции, протекающей в условиях, когда сжимаемость среды несущественна, т.е. для жидкостей, что собственно и положили Обербек и Буссинеск в основу своего вывода уравнений свободной конвекции из (1.1)–(1.3)
divv=0, (1.6) где р0 - плотность среды в невозмущенном состоянии; a, v = ju/р0,/3 коэффициенты температуропроводности, кинематической вязкости и объемного расширения жидкости; ё - единичный вектор, направленный противоположно вектору ускорения силы тяжести.
Существует и другой подход при рассмотрении свободной конвекции в условиях, когда размеры тела, нарушающего тепловое равновесие, малы по сравнению с объемом окружающей его жидкости. В этом случае область теплового и гидродинамического возмущения локализуется около рассматриваемого тела. Вне этого пограничного слоя жидкость можно считать неподвижной [4,5]. Наличие тонкого пограничного слоя позволяет существенно упростить систему уравнений (1.1)-(1.3), воспользовавшись стандартными приближениями пограничного слоя [6,7], которые обобщены на случай криволинейной смоченной границы в [8] путем кусочно-гладкой аппроксимации: где u,v - компоненты вектора скорости, t0 - характерная температура системы, (р - текущий угол между нормалью поверхности и вектором ускорения силы тяжести. При (р = 0 Польгаузеном [9] показана корректность такого подхода для описании свободноконвективного течения около полуограниченной вертикальной пластины при мгновенном изменении температуры стенки до постоянного значения, превышающего, например, температуру жидкости в объеме t0. Такой подход при моделировании свободноконвективных течений в замкнутых областях на примере вертикального цилиндра и сферы применен в [10,11]. Уравнения Навье-Стокса, составляющие основу уравнений свободной конвекции, обладают рядом специфических особенностей, одной из которых является пространственно-эллиптический характер решений [12], обусловленный влиянием вязкости во всем поле течения. В связи с этим для их решения требуется использовать типичные для эллиптических уравнений методы анализа. В отличие от уравнений пограничного слоя при этом требуется постановка граничных условий на всех границах рассматриваемой области. Кроме того, система уравнений Навье-Стокса нелинейна. Эта нелинейность, типичная для системы гидродинамического типа, обусловлена в случае несжимаемой жидкости инерционными составляющими в уравнениях количества движения. Уравнения конвекции в приближении Обербека-Буссинеска отличаются своей спецификой ввиду значительного взаимного влияния полей течения и температуры. В связи с этим нестационарность течения, обусловленная его неустойчивостью, обнаруживается для такого класса течений при меньших значениях гидродинамических критериев (например, числа Рейнольдса), чем в случае течения изотермической жидкости.
Достаточно полный анализ современного состояния проблемы адекватного моделирования в задачах свободной конвекции приведен в монографии [13]. Каждая из трех указанных моделей обладает преимуществами и недостатками, поэтому использование одной из них диктуется предметным приложением и геометрическими особенностями [14]. Так, общие уравнения Навье-Стокса наиболее точно описывают гидротермическую обстановку, но не обладают инвариантностью в отличие от уравнений Обербека-Буссинеска. Для каждой конкретной среды необходимо знание зависимостей теплофизических параметров от температуры, что, собственно говоря, требует заново формулировки и записи уравнений и поиска новых методов решения. Уравнения Обербека-Буссинеска, несмотря на универсальность, справедливы при достаточно небольших изменениях теплофизических параметров, т.е. фактически для неинтенсивных гидродинамических режимов и малых температурных градиентов. Использование погранслойного подхода теряет свою эффективность в замкнутых областях течения конечного размера.
Численное интегрирование
Среди различных способов получения новых знаний о закономерностях свободной конвекции главенствующую позицию занимает эксперимент. Также он служит критерием адекватности предлагаемых гипотез. Однако возникает ряд трудностей при его реализации, поскольку без внесения возмущений не удается снять измерения в исследуемых гидротермических полях, а бесконтактные методы чрезвычайно сложны и неоднозначны для математической интерпретации.
Задачей эксперимента являются измерения полей температуры и скорости. Обычно эти измерения классифицируют на локальные, полевые и интегральные [33]. В первом случае температура и скорость определяются в течение малого промежутка времени для малого объема жидкости. Как правило, фиксируется некоторая точка пространства и регистрируется зависимость скорости и температуры от времени. Локальные измерения выполняются с помощью датчиков, преобразующих измеряемые величины в электрические или оптические сигналы. Для определения полей скорости и температуры с помощью локальных датчиков необходимо сканирование исследуемого объема. Основной трудностью, которая возникает при реализации такой стратегии, является нивелирование возмущающих факторов от внесения датчиков в исследуемую область [34]. Если поля стационарны, то сканирование может выполняться в течение длительного времени, что не вызывает технических затруднений и позволяет оценить возмущающие факторы. Сканирование нестационарных полей в лучшем случае удается выполнить только вдоль некоторых линий.
В полевых методах измерения выполняются одновременно для всех точек некоторого объема жидкости. Очевидно, что полевые измерения могут быть только бесконтактными. Типичным примером полевых методов являются теневые и интерферометрические измерения распределений плотности жидкости или газа в двумерных потоках [35-39]. Мгновенное распределение плотности регистрируется в этом случае на одном снимке. Распространенным полевым методом измерения скорости является фотографирование светорассеивающих частиц.
Интенсивные экспериментальные 24,96-\/х исследования свободно-конвективного течения
В интегральных методах измеряются величины, которые являются функционалами от распределенных скоростей и температур. К типичным примерам таких величин можно отнести полный тепловой поток через некоторую поверхность, результирующую силу, действующую на обтекаемое потоком тело, среднюю температуру объема жидкости. С увеличением числа измеряемых в опыте интегральных величин конвективное движение описывается ими с возрастающей подробностью. Если эти величины образуют непрерывные множество, зависящее от того же числа переменных, что и исследуемое поле, то можно поставить задачу о восстановлении поля по его интегральным характеристикам.
Гидротермические поля по данным Шмидта и Бекмана около вертикальной пластины получения Польгаузеном классического анализа явлений переноса у вертикальной пластины, внезапно нагреваемой до конечной температуры, в неограниченной вязкой среде [1]. Проверке такого анализа, как отмечается в [4], были посвящены работы Шмидта и Бекмана, в которых анализировалось распределение скоростей и температур в воздухе около нагретой вертикальной плиты (рис. 1.6). Было показано, что при Pr 1 основная область тепловых и гидродинамических возмущений при свободной конвекции действительно сосредоточена в относительно узком слое жидкости около поверхности теплообмена, причем давление в каждом его горизонтальном сечении равно гидростатическому давлению в невозмущенной среде.
Можно констатировать, что в настоящее время достигнут гораздо более скромный уровень в исследовании свободно-конвективных течений, возникающих под действием сил плавучести в неравномерно нагретой среде, по сравнению с вынужденно-конвективными течениями. Оценка состояния уровня экспериментальных исследований подобных течений свидетельствует, в частности, о существенной недостаточности, а нередко и противоречивости полученных данных. Причины отмеченного, в определенном смысле, «хронического отставания» уровня экспериментальных исследований свободно-конвективных течений от аналогичного уровня исследований вынужденных течений связаны не только с особенностями данного течения (сравнительно небольшой уровень средних скоростей). Большие трудности возникают и при создании собственно экспериментальных установок, способных обеспечить высокостабильный свободно-конвективный поток в течение достаточно больших промежутков времени. По данным [40] в мире насчитывалось не более пяти установок (Япония, Франция, США), способных генерировать свободно-конвективные течения в широких диапазонах изменения числа Грасгофа.
Все сказанное относится и к экспериментальным исследованиям для внутренних задач свободной конвекции. Общепринятая физическая гипотеза, описывающая качественную картину явлений переноса, изложена в [41]. При прогреве неподвижной жидкости для постоянной плотности теплового потока на смоченных поверхностях температура стенки сначала растет, а в массе жидкости распространяется температурная волна, которая через определенное время достигает центра. После этого разность температур стенки и центра стремится к постоянному значению. Такой режим конвекции принято называть квазистационарным. Поля температур в центральной и придонной частях емкости разделяются на ядро с равномерной температурой по горизонтали и пограничный слой у стенки, в котором происходит основное изменение температуры во времени. Следовательно, дополнительным критерием при расчете теплообмена внутри емкости является число Fo. Процесс теплообмена в емкости сопровождается тепловой стратификацией в ядре потока, при которой верхние слои жидкости имеют температуру выше средней. Стратификация может быть значительной при неравномерном нагреве емкости, тепловыделении в жидкости и тепломассообмене на свободной поверхности. При длительном нахождении жидкости в условиях притока теплоты необходимо учитывать возрастание давления. Такая гипотеза основывается на классических экспериментальных работах [42-45].
Следующий этап экспериментальных исследований сосредотачивался на детализации картины явлений переноса. Например, в [46] рассматривалась свободная конвекция в крупнотоннажном вертикальном резервуаре при охлаждении вязкой жидкости. Найдено хорошее совпадение экспериментальных данных по теплообмену, полученных на основе профилей температур и скоростей по объему, с результатами предыдущего этапа изучения на модельных установках. В [47] установлено, что известные модели турбулентности, например, к-є модель, не в состоянии адекватно предсказать переход от ламинарной к турбулентной конвекции, т.е. определить критическое число Грасгофа. Среди последних работ на эту тематику можно выделить [48], в которой проведен широкий эксперимент с изменением материала стенки, различных граничных условий и геометрических характеристик, однако новых результатов, кардинально отличающихся от предыдущих работ, не было получено. Работы [49,50] практически дублируют [48], и хотя их целью было существенное уточнение картины течения и переноса теплоты, тем не менее, принципиально новых результатов, которые существенно бы прояснили ситуацию с формулирование гидротермической структуры, тоже не был получено.
Решение задачи в прямоугольнике
Скорость. Численное исследование можно провести очень быстро. Конструктор имеет возможность меньше, чем за день, просчитать сотни вариантов и выбрать оптимальную конструкцию, в то время как соответствующее экспериментальное исследование заняло бы очень много времени.
Полнота информации. Численное решение задачи дает подробную и полную информацию. С его помощью можно найти значения всех имеющихся переменных (таких, как скорость, давление, температура, концентрация, интенсивность турбулентности) во всей области решения. В отличие от эксперимента для расчета доступна практически вся исследуемая область и отсутствуют возмущения процесса, вносимые датчиками при экспериментальном исследовании. Очевидно, что ни в одном экспериментальном исследовании невозможно измерить распределения всех переменных во всей исследуемой области. Поэтому, даже если проводится экспериментальное исследование, большое значение для дополнения экспериментальной информации имеют результаты численного решения.
Возможность математического моделирования реальных условий. Численное уравнение можно получить для реальных условий исследуемого процесса, что далеко не всегда возможно при экспериментальном исследовании. Возможность моделирования идеальных условий. Если с помощью численного решения изучаются закономерности физического процесса, а не сложные инженерные задачи, можно сконцентрировать внимание на нескольких существенных параметрах этого процесса и исключить все несущественные явления. При этом можно моделировать многие идеализированные условия, например двухмерность, постоянство плотности, адиабатическую поверхность или бесконечно быструю реакцию. При экспериментальном исследовании даже с помощью довольно тщательного эксперимента не всегда можно достичь таких идеализированных условий.
Недостатки численного решения. Численное решение дает количественное выражение закономерностей, присущих математической модели. Напротив, с помощью экспериментального исследования наблюдается сама действительность. Таким образом, полезность расчета ограничена обоснованностью математической модели. Следует заметить, что результат численного решения зависит как от численного метода, так и от математической модели. Если используемая математическая модель не соответствует изучаемому явлению, то с помощью даже очень хорошей численной методики можно получить не нужные результаты.
Для обсуждения недостатков численного исследования полезно разбить все практические проблемы на две группы.
Группа А: проблемы, для которых математическая модель достаточно обоснована (например, теплопроводность, ламинарные течения, простые турбулентные пограничные слои).
Группа В: проблемы, для которых обоснованные математические модели пока не разработаны (например, сложные турбулентные течения, течения некоторых неньютоновских жидкостей, образование окиси азота при турбулентном горении, некоторые двухфазные течения).
Конечно, ответ на вопрос, в какую группу попадет данная проблема, будет зависеть от того, какую математическую модель мы считаем обоснованной. Недостатки для группы А. Можно утверждать, что для большинства проблем группы А численное решение имеет очень большие преимущества перед экспериментальным исследованием. Однако если цель исследования очень узкая (например, надо определить падение давления в сложных аппаратах), численное решение может быть не дешевле, чем эксперимент. Для задач, включающих сложную геометрию, сильные нелинейности, значительное изменение свойств жидкости и т. д., получение численного решения может оказаться трудным и чрезмерно дорогостоящим, а то и невозможным. Наконец, если математическая постановка задачи допускает более одного решения, трудно определить, соответствуют ли результаты расчета действительности.
Недостатки для группы. Б. Проблемы этой группы разделяют все недостатки проблем группы А. Дополнительно к ним здесь неясно, в какой мере численные результаты согласуются с действительностью. В этих случаях требуется экспериментальное обоснование результатов численного исследования.
Вычислительный эксперимент в общем виде включает следующие этапы: 1 – выбор физической модели; 2 – построение математической модели, уравнений и краевых условий; 3 – разработка численного метода (схемы); 4 – разработка программного комплекса для решения задачи и обработки результатов; 5 – отладка и тестирование программной реализации; 6 – проведение вычислительного эксперимента; 7 – практическое применение результатов, сравнение их с данными лабораторного или натурного эксперимента.
Следующий за построением математической модели этап реализации численного эксперимента направлен на замену непрерывной области решения дискретной. Так как в данном случае речь идет о математических моделях в виде дифференциальных уравнений, то их численное решение всегда начинается с аппроксимации дискретными соотношениями, что подразумевает способ дискретизации непрерывной области решения [80,81]. Существуют два основных способа дискретизации: замена непрерывной области некоторой совокупностью точек (основа метода конечных разностей) и разбиение области решения совокупностью геометрических элементов (основа метода конечных элементов).
Вычисление нормы оператора перехода оптимального итерационного шага
Обозначим ДФ = Д1Ф + 2Д2Ф + ВАЗФ . Заменяя в уравнении (3.4) частные производные конечно-разностными аналогами, а известную функцию ее сужением на узлы сетки и отбрасывая слагаемые более высокого порядка малости, получим уравнение
Докажем следующее свойство оператора B. Теорема 6. Оператор B положительно-определенный. Доказательство. Поскольку B = B1 + 2B2+B3 и согласно теореме 1, линейная комбинация положительно-определенных операторов с положительными коэффициентами положительно определена, достаточно показать лишь положительную определенность операторов B 1,B2,B3.
Вгг+а2Вгг-1-аЛ-1=0 Здесь =Д(), а _1=А(1,2/ ) - минор, получаемый из определителя А(к) вычеркиванием первой строки и второго столбца. Исключая Вп_1 и Вп_2 из первого уравнения получим характеристическое
Вычислив 23 1 12 постоянные, получим dQt[Bl(n)] = \ + n + п2+п3+ п4. Поскольку dettS»] 0 \fn є N ,Bx{n)удовлетворяет критерию Сильвестра. Матрица В3 симметрична и в другом ортогональном базисе при замене / на j совпадает с матрицей Вх, то есть Въ и Вх подобны, поэтому Въ положительно определен. Матрица В2 представляет собой произведение положительно определенных квазидиагональной и квазипятидиагональной матриц. В силу перестановочности этого произведения согласно теореме 5, оно положительно определено. Таким образом, оператор В положительно определенный.
Учитывая граничные условия (3.14) - (3.15), рассмотрим спектральную задачу для оператора В в квадрате со стороной / = 1 - 2h.
Значение параметров а и J3 выберем так, чтобы функция u(x,y) = sin(ax)sm({3y) удовлетворяла граничным условиям задачи. При х = 0 или у = 0 граничное условие выполняется автоматически. При х = I
Норма оператора перехода вычисляется по формуле IIS113 =max//J, где fik - собственные значения, тогда из (3.19) следует
Таким образом, мы получили не только выражение для нормы оператора перехода, но и значение оптимального итерационного шага. В таблице 3.2 приведены расчеты входных параметров задачи (3.16) и показана сходимость схемы при отклонении от оптимального шага.
Чтобы доказать устойчивость схемы (3.16) воспользуемся теоремой [110], адаптированной к рассматриваемой задаче.
Теорема 7. Для устойчивости схемы (3.17) достаточно, чтобы выполнялось условие для разрешающего оператора при этом верна априорная оценка Поскольку оператор перехода схемы (3.17) постоянный с нормой меньшей единицы, очевидно неравенство: Таким образом схема (3.16) устойчива и выполняется оценка Далее рассматривается стационарное разностное уравнение, полученное из уравнения (3.16) отбрасыванием разностного аналога производной по фиктивному времени, в операторной форме ) Пусть и решение уравнения (3.23). Поскольку схема (3.16) точно аппроксимирует уравнение (3.23) на его решении, разность zk = Ok-u удовлетворяет однородному уравнению Решая его относительно zk+l, получим zk+l = Szk. Используя рекурсию, получим zk+l=Skz0, а поскольку Ц І І для любой размерности сетки, схема (3.16) сходится к решению стационарной задачи. Действительно O -M = ZJ = 5 Z0 5 -Z0 5 Z0 0, ОО. Количество шагов к для обеспечения точности є приближения к стационарному решению zk \\ є \\ z0 \\ определяется формулой: k , ) в которой ln(1S) - скорость сходимости. Поскольку погрешность аппроксимации задана неравенством (3.22), то и точность приближения к стационарному решению e = O(h2). В таблице 3 приведены расчеты количества шагов для обеспечения точности є = h2 и максимального по модулю значения искомой функции Фkh. Чтобы показать сходимость решения задачи Коши (3.16) к решению краевой задачи (3.4)-(3.7) воспользуемся теоремой [112], модифицированной под формулировку нашей задачи. Теорема 8. Если: 1) решение краевой задачи (3.4)-(3.7) существует в некотором классе функций, 2) разностная задача (3.16) аппроксимирует краевую (3.4)-(3.7) на классе решения, 3) разностная задача (3.16) корректна, то: при h —»0 решение Фh разностного уравнения стремится к решению Ф дифференциального уравнения, т.е. Ф-ФhH - 0. Существование и единственность решения краевой задачи (3.4)-(3.7) доказаны. Уравнение (3.16) аппроксимирует (3.4) со вторым порядком (3.22). Начальное (3.5) и граничное условия (3.6) аппроксимируются точно, граничное условие (3.15) аппроксимирует (3.7) с первым порядком точности. В силу определения пространства сеточных функции Hh разностная задача устойчива по граничному условию (3.7), поэтому корректность задачи следует из ее устойчивости по начальным условиям и устойчивости по правой части, которые были доказаны в теореме 7. Таким образом, все условия теоремы 7 выполнены, следовательно, решение разностной задачи (3.16) стремится к решению краевой задачи (3.4)-(3.7) при h - 0. На квадрате [0,1] х [0,1] введем сетку аналогично стационарному случаю. Будем рассматривать сеточные функции и введем гильбертово пространство Н м Д7 с нормами, определенными так же, как и в стационарной задаче (3.17). Конечно-разностная аппроксимация бигармонического оператора в каждый момент времени задана выражениями (3.9) - (3.11). Разностный