Рефераты - Афоризмы - Словари
Русские, белорусские и английские сочинения
Русские и белорусские изложения
 

Разработка и исследование цифровой модели теплового потока при течении вязкой жидкости в канале с внешними нагревающимися элементами

Работа из раздела: «Математика»

/

/

МИНИСТЕРСТВО ОБРАЗОВАНИЯ И НАУКИ РОССИЙСКОЙ ФЕДЕРАЦИИ

ФИЛИАЛ ФЕДЕРАЛЬНОГО ГОСУДАРСТВЕННОГО АВТОНОМНОГО ОБРАЗОВАТЕЛЬНОГО УЧРЕЖДЕНИЯ ВЫСШЕГО ПРОФЕССИОНАЛЬНОГО ОБРАЗОВАНИЯ

«НАБЕРЕЖНОЧЕЛНИНСКИЙ ИНСТИТУТ КАЗАНСКОГО (ПРИВОЛЖСКОГО) ФЕДЕРАЛЬНОГО УНИВЕРСИТЕТА» В Г.НАБЕРЕЖНЫЕ ЧЕЛНЫ

ФАКУЛЬТЕТ ПРИКЛАДНОЙ МАТЕМАТИКИ И ИНФОРМАЦИОННЫХ ТЕХНОЛОГИЙ

КАФЕДРА ПРИКЛАДНОЙ МАТЕМАТИКИ И ИНФОРМАТИКИ

ВЫПУСКНАЯ КВАЛИФИКАЦИОННАЯ РАБОТА

РАЗРАБОТКА И ИСЛЕДОВАНИЕ ЦИФРОВОЙ МОДЕЛИ ТЕПЛОВОГО ПОТОКА ПРИ ТЕЧЕНИИ ВЯЗКОЙ ЖИДКОСТИ В КАНАЛЕ С ВНЕШНИМИ НАГРЕВАЮЩИМИСЯ ЭЛЕМЕНТАМИ

Специальность: 010501.65 - Прикладная математика и информатика

Набережные Челны - 2012 год

Введение

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

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

Первая форма теплообмена представляет собой перенос тепла посредством теплопроводности. Она характеризуется тем, что ее возникновение обусловлено наличием вещественной среды, и тем, что теплообмен совершается только между непосредственно соприкасающимися частицами тела. Этот процесс можно себе представить как распространение тепла от частицы к частице.

Второй формой теплообмена является перенос тепла посредством конвекции. Этот перенос осуществляется в тех случаях, когда частицы тела меняют свое местоположение в пространстве и при этом выполняют роль носителей тепла. Такой процесс происходит в движущихся жидкостях и газах и всегда сопровождается передачей тепла посредством теплопроводности от частицы к частице при условии, что вся масса текущей жидкости не обладает одинаковой температурой. До тех пор пока мы рассматриваем области, расположенные внутри потока, и, следовательно, не обращаемся к процессам, происходящим на твердых поверхностях, которые ограничивают поток, или на его свободной поверхности, мы можем обе формы переноса тепла охватить одним понятием -- теплопроводность в движущихся средах. Если же рассматриваются и твердые стенки то в общем случае между стенками и движущейся средой наблюдается теплообмен, который обусловлен тем, что частицы среды, соприкасающейся со стенкой, воспринимают от нее тепло и уносят его с собой. Теплообмен между средой и стенкой называют теплоотдачей.

Третья форма теплообмена -- перенос тепла посредством излучении. Эта форма характеризуется тем, что некоторая часть внутренней энергии тела преобразуется в энергию излучения и уже в такой форме передается через пространство. Встречая на своем пути другое тело, лучистая энергия полностью или частично вновь превращается в теплоту.

За последние десятилетия интерес к математическому моделированию сложных физических процессов и необходимость в нем в нем заметно возросли. Этому в значительной мере способствует прогресс в развитии компьютерной техники, численных методов решения всех типов задач математической физики и реализуемых на этой основе математических моделей. Любая современная наукоемкая технология, так или иначе, использует результаты вычислительных экспериментов. Сказанное в полной мере относится к исследованиям в области динамики жидкости и тепломассообмена.

Для численного решения практических задач, связанных с теплопереносом, течением жидкости и другими аналогичными явлениями, требуется, как правило, интегрирование системы нелинейных дифференциальных уравнений в частных производных по пространственным координатам и времени. Хотя существуют численные методы для получения такого решения, задача написания и использования общих вычислительных программ для всех практически важных процессов тепломассопереноса достаточно трудна.

Важным классом задач, подобных задачам теплопроводности, являются задачи о полностью развитых течениях и теплопереносе в каналах. Вдали от входного участка обычно существует область, в которой продольная скорость и поле температуры демонстрируют некоторое особое поведение. В этой области поля скорости и температуры описываются уравнениями, похожими на уравнение, использующееся для задач двумерной теплопроводности.

По трубопроводам транспортируют газы и различные жидкости. Встречаются трубопроводы различных сечений, но чаще круглые и полукруглые. Часто к стенкам канала прикрепляют ребра для увеличения площади поверхности, за счет чего улучшается теплообмен. Присутствие ребер приводит также и к увеличения перепада давления по длине канала при том же расходе через его сечение.

Жидкость в трубопроводе подогревают для уменьшения ее вязкости и увеличения скорости перемещения. Важно проанализировать поведение жидкости в трубопроводе, изменение температуры от теплопроводности и количество прикрепленных ребер.

В данном случае исследование актуально именно потому, что определенные аспекты темы изучены не в полной мере и проведенное исследование направлено на преодоление этого пробела.

Целью работы является разработка и исследование цифровой модели теплового потока при течении вязкой жидкости в канале с внешними нагревающимися элементами.

Перед нами стоит ряд задач:

– Создание математической модели распределения температурного поля в вязкой жидкости;

– Выбор метода реализации математической модели;

– Разработка алгоритмов реализации модели с помощью выбранных методов;

– Выбор программной среды для его реализации;

– Разработка программы по алгоритму с помощью выбранного программного обеспечения;

– Разработка цифровой модели изменения поля температуры в зависимости от: теплопроводности жидкости и металла, граничных условий;

– Проверка адекватности модели.

1. Постановка задачи

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

Теплообмен между нагретым кожухом и жидкостью может происходить только через поверхность трубы.

Для исходной области заданы граничные условия на внешней и внутренней границе кожуха. На нижней части кожуха задается постоянная температура T1 (идет нагрев кожуха), на внутренней стенке кожуха задана начальная температура T3. Задается температура жидкости T2, которая втекает в область нагрева (Рисунок 2).

Задается радиус трубы R2 и кожуха R1. Для обеспечения жесткости и увеличения площади теплообмена к внутренней стороне трубы присоединяют ребра жесткости, длиной половины R2.

Задается теплопроводность материала кожуха k и вязкость жидкости м.

Из-за симметрии канала расчетной областью будет половина кольца (Рисунок 3).

Рисунок 3 - Расчетная область канала

Требуется исследовать, как выглядит температурное поле жидкости при изменении:

– температуры нагрева кожуха;

– теплопроводность материала кожуха;

– вязкости жидкости;

– конфигурацию трубы (изменение количество ребер);

– температуры стенки кожуха и жидкости.

2. Основные понятия теории течения жидкости

2.1 Общие характеристики жидкости

Жидкость -- одно из агрегатных состояний вещества, материал, обладающий свойством текучести, т.е. сплошная легкодеформируемая среда, обладающая свойством течь. Основным свойством жидкости, отличающим её от других агрегатных состояний, является способность неограниченно менять форму под действием касательных механических напряжений, даже сколь угодно малых, практически сохраняя при этом объём [4].

К основным параметрам, характеризующим свойства жидкости, принадлежат:

Вязкость м

Жидкости (как и газы) характеризуются вязкостью. Вязкостью, или внутренним трением жидкости, называется сопротивление, которое жидкость оказывает перемещению одного ее слоя относительно другого под влиянием внешней силы.

Механизм внутреннего трения в жидкостях и газах заключается в том, что хаотически движущиеся молекулы переносят импульс из одного слоя в другой, что приводит к выравниванию скоростей -- это описывается введением силы трения.

Различают динамическую вязкость (единицы измерения: пуаз, 0,1Па·с) и кинематическую вязкость (единицы измерения: стокс, мІ/с, внесистемная единица -- градус Энглера). Кинематическая вязкость может быть получена как отношение динамической вязкости к плотности вещества и своим происхождением обязана классическим методам измерения вязкости, таким как измерение времени вытекания заданного объёма через калиброванное отверстие под действием силы тяжести [3].

Плотность с

Плотностью жидкости или газа называется количество массы, заключающейся в единице объема. Для определения плотности вещества нужно массу тела m разделить на его объем [3].

Размерность плотности -- кг/м3 (система единиц СИ).

Теплоемкость Ср

Теплоёмкость -- физическая величина, определяющая отношение бесконечно малого количества теплоты дQ, полученного телом, к соответствующему приращению его температуры дT.

Единицей СИ для удельной теплоёмкости является джоуль на килограмм-цельсий. Следовательно, удельную теплоёмкость можно рассматривать как теплоёмкость единицы массы вещества. На значение удельной теплоёмкости влияет температура вещества [3].

Теплопроводность k

Теплопроводность -- это молекулярный перенос теплоты между непосредственно соприкасающимися телами или частицами одного тела с различной температурой, при котором происходит обмен энергией движения структурных частиц (молекул, атомов, свободных электронов) [2].

Существует также понятие идеальной жидкости.

Идеальная жидкость -- в гидродинамике -- воображаемая (идеализированная) жидкость, в которой, в отличие от реальной жидкости, отсутствует вязкость. В идеальной жидкости отсутствует внутреннее трение, то есть нет касательных напряжений между двумя соседними слоями.

Моделью идеальной жидкости пользуются при теоретическом рассмотрении задач, в которых вязкость не является определяющим фактором и ею можно пренебречь. В частности, такая идеализация допустима во многих случаях течения, рассматриваемых гидроаэромеханикой, и даёт хорошее описание реальных течений жидкостей и газов на достаточном удалении от омываемых твёрдых поверхностей и поверхностей раздела с неподвижной средой. Математическое описание течений идеальных жидкостей позволяет найти теоретическое решение ряда задач о движении жидкостей и газов в каналах различной формы, при истечении струй и при обтекании тел [8].

/

/

2.2 Общие характеристики течения в каналах

Хоть каналы могут иметь любую геометрическую форму, для простоты рассмотрим прямой канал прямоугольного поперечного сечения, изображенный на рисунке 2.2. Течение в основном направлено вдоль оси z, (назовем ее продольной координатой, в отличие от поперечных координат х и у). Течение вызывается градиентом давления ?р/?z, который обычно отрицателен. Давление практически постоянно в поперечном сечении и изменяется вдоль оси z.

Так как изменения скорости и температуры по оси z малы по сравнению с изменениями их по оси х или у, то очень часто пренебрегают вязким напряжением, вызванным градиентом ?w/?z, и переносом тепла за счет градиента ?Т/?z. При этом обычно член ?р/?z в уравнении для продольной составляющей скорости полагают равным ?р/?z, где р соответствует среднему значению давления в поперечном сечении. Течение в канале может быть стационарным или нестационарным [13].

2.3 Ламинарное движение

Это движение, называют также потенциальным (безвихревым) движением.

При наблюдении за движением жидкости в трубах и каналах можно заметить, что в одном случае жидкость сохраняет определенный строй своих частиц, а в других - перемещаются бессистемно.

Ламинарным называется слоистое течение без перемешивания частиц жидкости и без пульсации скорости и давления. При ламинарном течении жидкости в прямой трубе постоянного сечения все линии тока направлены параллельно оси трубы, при этом отсутствуют поперечные перемещения частиц жидкости. Ламинарное течение - это плавное, упорядоченное, регулярное движение, когда отдельные струйки жидкости, не перемешиваясь, как бы скользят друг по другу. При таком потоке существует лишь молекулярное трение между соседними струйками [13-16].

Ламинарное течение жидкости наблюдается при небольших скоростях ее движения. Внешний слой жидкости, примыкающий к поверхности трубы, в которой она течет, из-за сил молекулярного сцепления прилипает к ней и остается неподвижным (рисунок 2.3). Скорости последующих слоев тем больше, чем больше их расстояние до поверхности трубы, и наибольшей скоростью обладает слой, движущийся вдоль оси трубы.

Бугорки и шероховатости обтекаются ламинарным потоком и не оказывают влияние на сопротивление трубы. При увеличении скорости движения толщина ламинарного слоя уменьшается.

При увеличении скорости потока траектории отдельных струек приобретают волнообразный характер, через некоторое время струйка исчезает, перемешиваясь с жидкостью. Характер течения при этом изменился. Траектории отдельных частиц приобретают хаотичный неустановившийся характер. Такое течение называется турбулентным - хаотичным, вихревым [16].

Переход от ламинарного к турбулентному режиму наблюдается при определенной скорости движения жидкости. Эта скорость называется критической хкр. Для характеристики течения жидкости используют так называемое число Рейнольдса - Re [1]. С ростом Re в некоторый момент происходит потеря устойчивости движения, струйки перемешиваются, в потоке образуются хаотически пульсирующие вихри. Такое неупорядоченное движение вязкой жидкости становится турбулентным. При турбулентном режиме обычно резко растёт сопротивление потока жидкости.

Значение этой скорости прямо пропорционально кинематической вязкости жидкости и обратно пропорционально диаметру трубы.

где н - кинематическая вязкость; k - безразмерный коэффициент; d - внутренний диаметр трубы.

Входящий в эту формулу безразмерный коэффициент k одинаков для всех жидкостей и газов, а также для любых диаметров труб. Этот коэффициент называется критическим числом Рейнольдса Reкр и определяется следующим образом:

Как показывает опыт, ламинарное движение возможно при сравнительно невысоких числах Рейнольдса. Для труб круглого сечения число Рейнольдса Reкр примерно равно 2300. Для нашего случая, когда канал в сечении кольцо и поверхность контакта между жидкостью и стенкой больше, то и больше локальных возмущающих факторов число Рейнольдса будет Reкр =1600 [1].

Критерий подобия числа Рейнольдса позволяет судить о режиме течения жидкости в трубе. При малых значениях числа Рейнольдса () наблюдается ламинарное течение, переход от ламинарного течения к турбулентному происходит в области , а при (для гладких труб) течение -- турбулентное. Если число Рейнольдса одинаково для разных жидкостей, то и режим течения таких жидкостей в трубах разных сечений одинаков.

Формула, по которой определяется число Рейнольдса в канале любой геометрии [16],

где - скорость, - гидравлический диаметр, Р - смоченный периметр,

- кинематическая вязкость, - динамическая вязкость.

Гидравлический диметр при конфигурации канала, когда сечение представляет собой кольцо [15]:

где Р - смоченный периметр, D и d - внешний и внутренний диаметр трубы.

2.4 Теплообмен

Теплообмен -- это самопроизвольный необратимый перенос теплоты (точнее, энергии в форме теплоты) между телами или участками внутри тела с различной температурой. В соответствии со вторым началом термодинамики теплота переносится в направлении меньшего значения температуры. Теплообмен существен во многих процессах нагревания, охлаждения, конденсации, кипения и оказывает значительное влияние на массообменные процессы [19].

Движущиеся среды, участвующие в теплообмене и интенсифицирующие его, называются теплоносителями (обычно капельные жидкости, газы и пары). Известны два основных способа проведения тепловых процессов: путем теплоотдачи и теплопередачей.

Теплоотдача -- теплообмен между поверхностью раздела фаз (чаще твердой поверхностью) и теплоносителем.

Теплопередача -- теплообмен между двумя теплоносителями или иными средами через разделяющую их твердую стенку либо.

Различают три разных механизма распространения теплоты: теплопроводность, конвективный и лучистый перенос.

Теплопроводность -- перенос энергии от более нагретых участков тела к менее нагретым в результате теплового движения и взаимодействия микрочастиц (атомов, молекул, ионов и др.). В чистом виде теплопроводность может встречаться в твердых телах, не имеющих внутренних пор, и в неподвижных слоях жидкостей, газов или паров. Количество переносимой теплопроводностью энергии, определяемое как плотностью теплового потока qT [Вт/(м2* К)], пропорционально градиенту температуры (закон Фурье):

где л - коэффициент теплопроводности, характеризующий его способность проводить теплоту, Вт/(м*К); знак минус указывает направление переноса теплоты в сторону снижения температуры.

Численная характеристика теплопроводности материала равна количеству теплоты, проходящей через материал площадью 1 кв.м за единицу времени (секунду) при единичном температурном градиенте.

С повышением температуры теплопроводность жидкостей, за исключением воды, уменьшается, а для всех других тел увеличивается.

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

где w - вектор скорости текучей среды; с, С, Т - плотность, теплоемкость и температура среды.

Плотностью жидкости или газа называется количество массы, заключающейся в единице объема. Для определения плотности вещества нужно массу тела m разделить на его объем. Размерность плотности - кг/м3 (система единиц СИ).

Теплоёмкость -- физическая величина, определяющая отношение бесконечно малого количества теплоты дQ, полученного телом, к соответствующему приращению его температуры дT.

Единицей СИ для удельной теплоёмкости является джоуль на килограмм-цельсий. Следовательно, удельную теплоёмкость можно рассматривать как теплоёмкость единицы массы вещества. На значение удельной теплоёмкости влияет температура вещества.

В большинстве случаев значения w, с, С и Т потоков теплоносителей таковы, что в направлении движения конвективный перенос преобладает над теплопроводностью. Однако при малых скоростях течения высокотеплопроводных жидкостей (расплавов металлов) может наблюдаться обратное соотношение.

По мере приближения к твердой поверхности, где скорость вязких жидкостей стремится к нулю, qт и qk также становятся сравнимы по величинам. При ламинарном режиме течения в направлении, поперечном движению, конвективный перенос отсутствует.

3. Математическая постановка задачи

3.1 Обобщенная математическая постановка задачи

Применение первого закона термодинамики к бесконечно малому контрольному объему в твердом теле или в неподвижной среде приводит к известному уравнению теплопроводности [1].

Обозначив три пространственные переменные через декартовы координаты х, у, z и время через t, можем записать уравнение теплопроводности в виде

где с -- плотность, с -- удельная теплоемкость.

В этом уравнении величины с, с, k могут зависеть от х, у, z, t и быть функциями температуры Т.

Уравнение теплопроводности является прямым следствием закона сохранения, представленного первым законом термодинамики, и пропорциональности плотности потока градиенту температуры. Существует множество других физических процессов, при которых соответствующая плотность потока некоторой величины пропорциональна градиенту этой величины и для которых существует закон сохранения. Отсюда следует, что эти процессы будут описываться дифференциальными уравнениями, аналогичными (3.1).

К подобным процессам можно отнести диффузию химических компонент, движение заряженных частиц в электромагнитном поле, течение в пористых материалах, потенциальные течения, перенос тепла и влаги в почве, а также полностью развитые течение и теплообмен в каналах. Построив вычислительную процедуру для решения уравнения (3.1), мы сможем применить ее и для любого аналогичного процесса, просто придавая новый смысл величинам Т, k и др. Например, можно интерпретировать Т как концентрацию, k как коэффициент диффузии и т.п. Удобнее работать с таким обобщенным дифференциальным уравнением, так как уравнение теплопроводности и другие аналогичные уравнения станут его частными случаями. В дальнейшем будем основываться на подобном обобщенном дифференциальном уравнении [2].

3.2 Дифференциальное уравнение для обобщенной переменной

Обозначим через ц зависимую переменную в обобщенном дифференциальном уравнении. Градиент ц приводит к соответствующему диффузионному потоку. Таким образом, плотность диффузионного потока Jx в направлении х задается в виде:

где Г -- обобщенный коэффициент диффузии.

Обобщенное дифференциальное уравнение, представляющее закон сохранения для ц, может быть записано в виде

где л -- величина, подобная объемной теплоемкости.

Сравнивая уравнение теплопроводности (2.1) с обобщенным дифференциальным уравнением (2.3), заключаем, что для того чтобы обобщенное уравнение описывало процесс теплопроводности, нужно сделать следующую замену: ц = Т, л = сс, Г = k,

В общем случае для значения ц существуют соответствующие значения величин л, Г. В простой задаче эти величины могут быть постоянными, но в общем случае нет ничего необычного в том, что они зависят от пространственных координат и времени, а также от самой переменной ц.

Обобщенное дифференциальное уравнение (2.3) записано в декартовой системе координат. Более компактно оно может быть записано:

где применено соглашение о суммировании по повторяющемуся индексу. Например, для трехмерного случая можно записать

здесь х1, х2 и х3 -- декартовы координаты (х, у и z) по трем направлениям. Будем называть три члена уравнения (2.5) нестационарным, диффузионным и источниковым членом соответственно.

В более общем виде выражение для плотности диффузионного потока (2.2) может быть записано следующим образом [3]:

где Ji -- плотность диффузионного потока в направлении координаты хi.

Диффузионный член в (2.5) удобно представить в виде:

Для двумерной задачи, сформулированной в полярных системах координат, выражение для будет иметь вид:

Плотности потоков соотносятся с градиентом ц следующим образом:

Уравнение для стационарного медленного течения в канале записывается в виде

где ср -- удельная теплоемкость при постоянном давлении; k -- теплопроводность жидкости. Левая часть (2.10) соответствует конвективному переносу в канале, а члены правой части описывают теплопроводность в жидкости. Так как обычно перенос тепла вдоль оси z очень мал по сравнению с переносом в поперечном его сечении, то последним членом в уравнении (2.10) можно пренебречь.

При простом полностью развитом течении в канале поперечные скорости и и v равны нулю, поэтому (2.10) упрощается и принимает вид

Уравнение (2.11) соответствует стационарной форме обобщенного дифференциального уравнения при следующей замене:

4. Методы решения

4.1 Выбор метода решения

Известны [16] несколько основных типов численных методов, применяемых для решения уравнений теплопроводности:

– Метод конечных разностей;

– Метод контрольного объема;

– Метод конечных элементов.

Получить аналитическое решение уравнений описывающих явления гидродинамики, как правило, не представляется возможным. В связи с этим используются численные методы. Все методы включают два этапа.

– На первом этапе все дифференциальные уравнения сводятся к алгебраическим уравнениям относительно значений искомых переменных (скорости, давления, температуры и т.п.) в конечном числе точек в пределах области решения. Этот этап называют дискретизацией уравнений

– На втором этапе производится решение алгебраических выражений с помощью соответствующего численного метода.

Различные типы численных методов вычислительной гидродинамики используют свои способы дискретизации основных уравнений.

Каждый из упоминавшихся выше численных методов имеет свои преимущества и недостатки. Преимуществами метода контрольных объемов являются:

· Простота применения к расчетным областям произвольной формы

· Наличие отработанных алгоритмов

· Ошибка дискретизации падает с ростом числа вершин

· Основное внимание уделяется балансу потоков через контрольные объемы

· Автоматическое выполнение законов сохранения

Граничные условия

Для численного моделирования требуется задание значений зависимых переменных или их нормальных градиентов на границах расчетной области. Граничные условия для нашей задачи были определены при постановке задачи.

Решение дискретных аналогов уравнений не может быть получено, если граничные условия не заданы. Граничные условия делятся на два типа - с известными значения на границе (задача Дирихле) и известными нормальными градиентами на границе (задача Неймана).

Начальные условия

Для нестационарных задач перед началом расчета необходимо задать значения всех искомых переменных во всей области потока в начальный момент времени.

Источники ошибок

В общем случае существует три основных первопричины ошибок:

· Допущения, принимаемые при построении математической модели

· Аппроксимации при дискретизации уравнений

· Отсутствие сходимости в процессе решения

О принимаемых допущениях говорилось ранее - они связаны с неточным описанием реальных физических явлений используемыми уравнениями. К этой же категории можно отнести аппроксимацию геометрии и граничных и начальных условий.

Ошибки дискретизации являются результатом применения аппроксимирующих соотношений для расчетов потоков через грани ячеек и источниковых членов в пределах ячеек.

Под отсутствием сходимости в данном случае подразумевается невозможность получения решения дискретных алгебраических уравнений с желаемой точностью. (Стремление численного решения к точному при последовательном измельчении сетки также часто называется сходимостью).

жидкость математический температура цифровой

4.2 Методы решения дифференциальных уравнений

Для решения задач математической физики, не поддающихся аналитическим решениям, используют численные методы, позволяющие аппроксимировать исходные дифференциальные уравнения. Развитие ЭВМ привело к появлению различных численных методов решения уравнений математической физики. Кроме того, за время развития вычислительных методов наметилась их специализация по отраслям приложения. Так, например, при моделировании глобальных атмосферных процессов и при решении акустических задач отдается предпочтение спектральным методам; методы конечных элементов применяются сейчас, в основном, для решения задач механики сплошной среды. Для решения проблем теплообмена и механики жидкости и газа используются как метод конечных разностей, так и метод контрольных объемов [12].

4.2.1 Метод конечных разностей

На рисунке 4.1 изображено разбиение одномерной расчетной области на узловые точки.

/

/

Конечно-разностный метод основан на замене производных дифференциального уравнения физического процесса в точке на конечные разности. Производные аппроксимируются разложением в ряды Тейлора. Для аппроксимации производной в точке 2 разложим исследуемую переменную Ц в ряд Тейлора в точках 1 и 3

Пренебрегая остаточным членом Rn(x) обоих рядов, вычитая и складывая уравнения, получим:

Заменяя производные в исходном дифференциальном уравнении в соответствии с (2.6), получаем разностный аналог, записываемый в виде системы алгебраических уравнений.

Данный метод не учитывает физической картины процесса, а наличие остаточного члена в ряде Тейлора и предположение о полиномиальном характере изменения зависимой переменной, приводит к нежелательным последствиям, например, для случая экспоненциального изменения переменной. Этим методом решается значительное количество практических задач газовой динамики и диффузии (теплообмена), однако в ряде случаев конечно-разностный метод не подходит для решения уравнений математической физики. Так, для решения дифференциальных уравнений с разрывными коэффициентами в конечно-разностном методе применяют схему с искусственной вязкостью, метод характеристик и дивергентные схемы. Указанные разновидности конечно-разностного метода применяются для задач гидро- и газодинамики [2].

4.2.2 Метод сеток

Дифференциальные уравнения вида решаются преобразованием их в алгебраические уравнения, называемые дискретными аналогами. Эти уравнения содержат в качестве неизвестных значения ц в выбранных дискретных местоположениях, которые образуют сетку и называются расчетными точками. Вокруг каждой расчетной точки строится контрольный объем, и дискретные аналоги получаются интегрированием уравнения по таким контрольным объемам.

Основная идея метода контрольного объема легка и понятна и поддается прямой физической интерпретации. Расчетную область разбивают на некоторое число контрольных объемов, грани которых показаны штриховыми линиями. После в геометрический центр каждого контрольного объема помещается узловая точка.

Видно, что некоторая расчетная точка «сообщается» с четырьмя соседними через четыре грани контрольного объема. Одна из граней приграничного контрольного объема совпадает с границей расчетной области, а граничная точка помещена в центр грани контрольного объема.

Удобно представлять контрольный объем нулевой толщины для граничной точки.

Границы контрольных объемов должны располагаться так чтобы, они совпадали с разрывами в свойствах материала, источниковых членах, граничных условиях и др.

Существует большое число геометрических величин, которые имеют отношение не к расчетным точкам, а к граням контрольных объемов. Грань с номером i лежит между точками i-1 и i. Другими словами, грань имеет тот же номер, что и ближайшая к ней точка в положительном направлении оси координат, т.е. в направлении увеличения i или j. Частным следствием такого построения контрольных объемов и используемой нумерации является то, что грань I = 2 и точка I = 1 совпадают с левой границей расчетной области.

Границы контрольных объемов должны располагаться так чтобы, они совпадали с разрывами в свойствах материала, источниковых членах, граничных условиях и др.

В соответствии с местоположение расчетных точек (рисунок 4.2) в направлении осей параметры и и у определяются индексами I и J. Индекс I увеличивается вдоль оси и, а индекс J -- вдоль оси у. Значение I=1 задает сеточную линию на левой границе, в то время как I=L1 указывает на правую границу. Аналогично J=1 и M1 соответствуют нижней и верхней границам. На самом деле во многих местах в программе левая и нижняя граница области называются границами I1 и J1 соответственно, а правая и верхняя -- границами L1 и M1. Логичность такого обозначения очевидна. Для удобства введены некоторые дополнительные переменные, определяемые так:

Таким образом, сеточные линии I=2 и L2 являются первыми внутренними линиями сетки вблизи соответствующих границ, а линии I=3 и L3 -- вторыми внутренними линиями. Подобными свойствами обладают линии J=M2 и МЗ. Таким образом, в диапазонах I=2, L2 и J=2, М2 задаются все внутренние расчетные точки. Зависимая переменная, например температура, в расчетной точке [I, J] будет обозначаться как Т[I, J]. Граничные значения Т обозначены как Т[1, J], Т[L1, J], Т[I, 1] и Т[I, M1].

Программа разработана для полярной системы координат и - продольная координата, y,r - координаты в радиальном направлении. Разница между ними заключается в том, что координата r должна измеряться от полюса, в то время как координата y может быть выбрана произвольно.

Существует большое число геометрических величин, которые имеют отношение не к расчетным точкам, а к граням контрольных объемов. Грань с номером i лежит между точками i-1 и i. Другими словами, грань имеет тот же номер, что и ближайшая к ней точка в положительном направлении оси координат, т.е. в направлении увеличения i или j. Частным следствием такого построения контрольных объемов и используемой нумерации является то, что грань I = 2 и точка I = 1 совпадают с левой границей расчетной области.

Дифференциальные уравнения законов сохранения массы, энергии, количества движения интегрируют по каждому объему. В результате находят дискретный аналог дифференциального уравнения, в который входят значения в нескольких узловых точках.

Полученный подобным образом дискретный аналог выражает закон сохранения для конечного контрольного объема точно так же, как дифференциальное уравнение выражает закон сохранения для бесконечно малого контрольного объема. Одним из важных свойств метода контрольного объема является то, что в нем заложено точное интегральное сохранение таких величин, как масса, количество движения и энергия на любой группе контрольных объемов и, следовательно, на всей расчетной области. Это свойство проявляется при любом числе узловых точек, а не только в предельном случае очень большого их числа [1].

Именно этот метод и будет использован в дальнейшей работе.

4.3 Полярная система координат

Поставленная задача решена с применением полярной системы координат. Все вычисления произведены в полярных координатах без перевода их в декартовые.

Полярная система координат -- двумерная система координат, в которой каждая точка на плоскости определяется двумя числами -- полярным углом и полярным радиусом. Полярная система координат особенно полезна в случаях, когда отношения между точками проще изобразить в виде радиусов и углов; в более распространённой, декартовой или прямоугольной системе координат, такие отношения можно установить только путём применения тригонометрических уравнений.

Полярная система координат задаётся лучом, который называют нулевым или полярной осью. Точка, из которой выходит этот луч, называется началом координат или полюсом. Любая точка на плоскости определяется двумя полярными координатами: радиальной и угловой. Радиальная координата (обычно обозначается r) соответствует расстоянию от точки до начала координат. Угловая координата, также называется полярным углом или азимутом и обозначается ц, равна углу, на который нужно повернуть против часовой стрелки полярную ось для того, чтобы попасть в эту точку.

Определённая таким образом радиальная координата может принимать значения от нуля до бесконечности, а угловая координата изменяется в пределах от 0° до 360°. Однако, для удобства область значений полярной координаты можно расширить за пределы полного угла, а также разрешить ей принимать отрицательные значения, что отвечает повороту полярной оси по часовой стрелке.

Одной из важных особенностей полярной системы координат является то, что одна и та же точка может быть представлена бесконечным количеством способов. Это происходит потому, что для определения азимута точки нужно повернуть полярную ось так, чтобы он указывал на точку. Но направление на точку не изменится, если осуществить произвольное число дополнительных полных оборотов.

Углы в полярных координатах задаются либо в градусах, либо в радианах. Выбор, как правило, зависит от области применения. В навигации традиционно используют градусы, в то время как в некоторых разделах физики, и почти во всех разделах математики используют радианы.

4.4 Обобщенное дискретное уравнение

Значения обобщенной зависимой переменной ц сохраняются в расчетных точках. Дискретное уравнение связывает значение ц в одной расчетной точке со значениями ее в четырех соседних точках. Такое уравнение получается при интегрировании уравнения по контрольному объему, содержащему точку (i, j).

Типичный контрольный объем показан на рисунок 4.5. Проинтегрировав по нему уравнение (2.2), получим

В формуле (4.1) верхний индекс «0» обозначает известное значение ц в начале шага по времени Дt; J -- плотность диффузионного потока через грань контрольного объема; -- осредненный по контрольному объему источниковый член; A, ДV -- площадь грани и контрольный объем.

Диффузионные потоки па гранях контрольного объема e и w могут быть рассчитаны следующим образом:

где De -- проводимость между точками P и E, которая вычисляется по значениям Г в этих точках. Зная размеры дxе- и дxе+, можно рассчитать De, по формуле

Аналогичным образом определяются значения D на других гранях. Проводимость зависит от площади грани контрольного объёма .

Так как источниковый член S может зависеть от ц, желательно учесть эту зависимость в линеаризованной форме. С этой целью S записывается в виде

где SP -- коэффициент при цP; Sc -- постоянная часть, которая от цР явно не зависит.

Когда нет необходимости в линеаризации, следует положить Sp равным пулю, а Sc приравнять к . В любом случае Sp никогда не должно иметь положительного значения.

Подстановка вышеприведенных выражений для J и S в (4.1) приводит к окончательной форме дискретного аналога. Он записывается в виде

где

Здесь проводимости De, Dw, Dn и Ds определяются в виде (4.3).

4.5 Представление граничных условий

Для каждой приграничной расчетной точки соответствующая граничная точка является соседней. Поэтому для граничной точки должно быть задано или значение зависимой переменной, или уравнение для его отыскания. Так как для всех четырех границ расчетной области граничные условия реализуются одинаковым образом, рассмотрим только левую границу. Представим два типа реализации граничных условий.

4.5.1 Первый порядок аппроксимации

На рисунке 4.6 показана окрестность граничной точки (i,j). Грань контрольного объема при i = 2, совпадающая с левой границей области, может рассматриваться лежащей между точками с ц1 и ц2, если вокруг точки с ц1 представить контрольный объем бесконечно малой толщины. Тогда плотность потока J2 на левой границе задается формулой (4.2). Удобно переписать ее в виде

На рисунке 4.6 видно, что описанный выше способ вычисления плотности потока J2 на границе соответствует односторонней схеме, так как грань лежит не посередине между точками с переменными ц1 и ц2, что дает не очень точные результаты. Так как представление граничных условий сильно влияет на все решение и значения плотностей потоков на границе часто являются важным результатом расчетов, то желательно получить более точную формулу для их определения.

4.5.2 Второй порядок аппроксимации

Формула более высокого порядка аппроксимации может быть получена, если считать, что плотность диффузионного потока меняется линейно между гранями контрольного объема. Предполагаемое распределение J между точками 1 и 2 описывается формулой

где

Проинтегрируем выражение (4.14) по х для получения профиля ц. Интегрирование производится только по области, лежащей между точками 1 и 2, где Г имеет постоянное значение. Возникающая при интегрировании константа находится из условия

ц = ц1 при х - х1 = 0

ц = ц 2 при х - х1, = д.

Использование этого условия после некоторых алгебраических преобразований приводит к выражению

Обобщенная форма этого выражения (4.18) имеет вид

Если множитель в = 4/3, получаем выражение (4.18). При в = 1 это выражение приводит к аппроксимации первого порядка точности, заданной формулами.

Выражение для J3 с учетом (4.2) может быть записано в виде

где D3 -- соответствующая проводимость между точками с переменными ц 2 и ц3.

Дискретный аналог для приграничного контрольного объема с индексами узловой точки (2,J) должен быть основан на выражении (4.19) для плотности потока J2 через его левую грань, а не на (4.19). В результате коэффициенты aW(2) и aE(2) записываются в виде:

Исходя из (3.12), выражение для плотности диффузионного потока J2 может быть записано в виде

где

Итак, выражения (4.23) -- (4.25) могут быть использованы для обоих вариантов аппроксимации граничных условий при соответствующем выборе значении в:

второй порядок аппроксимации при в = 4/3;

первый порядок аппроксимации при в =1.

Индикатор граничных условий.

Тип граничных условий на левой, правой, нижней и верхней границах расчетной области определяется индикаторами KBC. В случае КВС=1 на границе известно значение ц.

Когда известно значение зависимой переменной цв на границе, то необходимо приравнять ее к постоянному члену b дискретного уравнения. Таким образом, если KBC=1, то необходимые изменения заключаются в следующем:

5. Компьютерные методы решения

5.1 Выбор и обоснование программного инструмента

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

Такие требования к современному программированию привели к созданию многочисленных RAD-систем (от англ. RAD -- Rapid Application Development -- быстрая разработка приложений), представляющих собой интегрированные среды разработчика, включающие в себя:

§ средства быстрого и удобного построения программ, в том числе визуального;

§ встроенные компиляторы и отладчики;

§ системы коллективной разработки проектов и т.д.

Одной из таких RAD-систем является Delphi. Итак, Delphi -- это объектно-ориентированная среда для визуального проектирования Windows приложений с развитыми механизмами повторного использования программного кода. Основным конкурентом Delphi является среда разработки Microsoft Visual C++, имеющая свои преимущества и недостатки, однако являющаяся более популярной, в основном, в силу того, что разработана именно фирмой Microsoft. Существенной чертой Delphi является компонентная модель разработки программных продуктов. Суть модели заключается в поддержке системой постоянно расширяемого набора объектных компонентов, из которых и строится программа. Компоненты в Delphi просты для использования и развития, как результат сокрытия значительной части той структуры программы, которая близка к взаимодействию с операционной системой. Таким образом, для создания в Delphi несложных программных продуктов совершенно не обязательно понимать внутреннюю структуру Windows-приложения, получаемого после разработки в Delphi. Достаточно просто уметь работать с некоторыми компонентами, поставляемыми вместе со средой разработчика. При этом начать работу со средой можно практически без предварительного ознакомления, а написание первого приложения не потребует углубления в особенности системы. Этому отчасти способствует удобный интерфейс среды разработчика, не перегруженный излишними вопросами к разработчику [13].

Однако такой подход совершенно неприемлем для серьезного программирования, и, рано или поздно, придется освоить и основы программирования под ОС Windows, и серьезно изучить саму среду разработки Delphi, а также возможности, которые она предоставляет. Кроме того, конечно же, для создания качественных программных продуктов необходимо глубокое понимание компонентной модели.

Для разработки нашего проекта была выбрана именно среда программирования Borland Delphi. Среда визуального программирования, как уже упомянули выше, выбрана потому, что она обеспечивает возможность максимально простого, быстрого и качественного создания стандартных фрагментов программы.

5.2 Разработка алгоритма

Схема решения задачи:

1. Задание значений переменных

2. Предоставление информации о расчетной сетке

2.1. Расчет координат граней контрольных объемов

2.2. Расчет координат расчетных точек

xi:=0.5*(xui+1+xui)

yj:=0.5*(yvj+1+yvj)

3. Задание теплового потока qw и параметров жидкости:

с плотность (кг/м3),

cp теплоемкость (дж/кг*К),

k теплопроводность (вт/К*м),

м вязкость (Па*с).

4. Расчет коэффициентов дискретных аналогов

5. Решение системы алгебраических уравнений методом переменных направлений

6. Вычисление средней температуры и вывод значений

Если Тср<Тзад

то

7. Вывод таблиц распределения температуры и скорости

Иначе идем на п.4.

5.3 Пример реализации алгоритма

Пример расчета координат точек:

procedure Zgrid; // Область делится по осям x и y на различные зоны и для каждой зоны создается сетка

var

// задаем переменные

begin

// Построение сетки зона за зоной

// Рассматривается направление X

xu[2]:=0.; // координаты граней по х

i2:=2;

for nz:=1 to nzx do // по зонам

begin

fcvlx:=ncvx[nz];

ilast:=i2;

i1:=ilast+1;

i2:=ilast+ncvx[nz];

for i:=i1 to i2 do

begin

dd:=(i-ilast)/fcvlx;

if (powrx[nz]>0) then

xu[i]:=xu[ilast]+xzone[nz]*power(dd,powrx[nz])

else

xu[i]:=xu[ilast]+xzone[nz]*(1.-power((1.-dd),(-powrx[nz])));

end;

end;

L1:=i2;

// Рассматривается направление Y

yv[2]:=0.;

j2:=2;

for nz:=1 to nzy do

begin

fcvly:=ncvy[nz];

jlast :=j2;

j1:=jlast+1;

j2:=jlast+ncvy[nz];

for j:=j1 to j2 do

begin

dd:=(j-jlast)/fcvly;

if (powry[nz]>0) then

yv[j]:=yv[jlast]+yzone[nz]*power(dd,powry[nz])

else

yv[j]:=yv[jlast]+yzone[nz]*(1.-power((1.-dd),(-powry[nz])));

end;

end;

M1:=j2;

end;

procedure Phi; // вывод на печать информацию о сетке и двумерные поля зависимых переменных

var

i,j:integer;

begin

begin

for j:=2 to M2 do

for i:=2 to L2 do

begin

gam[i,j]:= cond*f[i,j,2]; // вычисление коэффициента диффузии

sc[i,j]:=-dpdz;

end;

end;

// определяем теплоизолированные части

for j:=2 to M2 do // задаем теплоизоляцию кроме первого сектора

begin

kbci1[j]:=2;

kbcl1[j]:=2;

end;

for i:=1 to ncvlx do

begin

kbcm1[i]:=2;

end;

end;

5.4 Разработка программы по алгоритму с помощью выбранного программного обеспечения

5.4.1 Описание некоторых процедур

Процедура READY

Функция READY заключается в расчете и сохранении значений геометрических величин. В памяти сохраняются только одномерные геометрические характеристики. Величины вида YCVR(J) *XCV(I), хотя и полезны для дальнейшей работы, не сохраняются, так как для этого требуется введение двумерного массива. Унифицированное использование двух систем координат (декартовой и полярной) становится возможным благодаря должному определению обобщенных геометрических характеристик, сделанному в READY. Ни в одном месте программы не нужно использовать оператор IF при программировании различных формул для каждой системы координат.

READY выполняет также еще одну функцию. Эта подпрограмма печатает заголовок с указанием используемой системы координат. Этот заголовок довольно важен, так как, если вы случайно используете не ту систему координат, сразу будете предупреждены. READY распечатывает текстовую переменную HEADER, которая содержит представленный пользователем заголовок задачи длиной до 64 символов.

Процедура HEART

Самая важная вычислительная часть программы содержится в подпрограмме HEART. Здесь рассчитываются коэффициенты дискретных аналогов и производятся модификации граничных условий. Взгляд, брошенный на подпрограмму HEART, прежде всего наткнется на внешний оператор цикла DO с параметрами N=1, NFMAX, предназначенный для организации последовательных вычислений каждой зависимой переменной. Внутри цикла производится вызов подпрограммы PHI, где пользователь определяет значения ALAM (I, J), GAM(I, J), SC (I, J), SP (I, J) и индикаторы граничных условий КВС. В заключение вызывается подпрограмма SOLVE для итерационного решения системы алгебраических уравнений.

В конце HEART увеличивается значение ITER, а если ITER достигает своего предельного значения, определяемого переменной LAST, то KSTOP присваивается значение единицы.

Процедура SOLVE

Задача подпрограммы SOLVE заключается в итерационном решении системы линейных алгебраических уравнений вида (2.11). Алгоритм решения, реализованный в SOLVE, описан ранее. Кроме повторов алгоритма нужное число раз в конце SOLVE реализованы некоторые дополнительные функции. Так, например, производится вычисление неизвестных значений ф на границе, когда индикатор граничных условий КВС = 2. Следует отметить, что алгоритм решения в SOLVE оперирует значениями ф только во внутренних точках. Когда решение получено, неизвестные значения ф на границе находятся из уравнений вида (2.40).

Как уже известно, в подпрограмме DEFLT большому числу переменных по умолчанию присваиваются некоторые значения. Однако некоторые переменные должны после каждого их использования принимать первоначальное значение. Такими переменными являются, например, источниковые члены SC(I, J) и SP(I, J), все параметры граничных условий КВС, и др. Этим переменным присваиваются их значения по умолчанию в конце подпрограмме SOLVE.

Процедура ZGRID

Область делится по осям х и у на различные зоны и для каждой зоны создается сетка.

Для нашей вычислительной процедуры очень важно, чтобы разрывы в распределении теплопроводности, источниковых членов и в граничных условиях совпадали с гранями контрольных объемов. При произвольном расположении разрывов не всегда можно добиться этого при использовании равномерной сетки или сетки, рассчитываемой по (3.1) и (3.2). В этом случае можно разделить расчетную область по оси х (так же, как и по оси у) на различные зоны таким образом, чтобы их границы совпадали с разрывами. Тогда можно задавать число контрольных объемов и значение п для каждой зоны в отдельности. Процедура ZGRID обеспечивает построение именно такой сетки.

Для ZGRID требуются следующие входные данные: число зон NZX и NZY по осям х и у соответственно, а также размер зоны, число контрольных объемов и значения п для каждой зоны в отдельности. Используя эти данные, ZGRID рассчитывает XU (I), YV (J), LI и Ml.

Пользователь должен задать значения массивов, характеризующих положения граней контрольных объемов XU (I) и YV (J). В общем случае сетка может быть неравномерной и конкретное распределение величин XU (I) и YV (J) зависит от особенностей задачи. Однако часто бывает необходимо задать равномерную сетку так, чтобы все контрольные объемы имели одинаковые размеры как по оси х, так и по оси у.

Входные данные, необходимые для ZGRID, состоят из значений переменных XL и YL, соответствующих размерам расчетной области по осям х и у, и из желаемого числа контрольных объемов. При этом с помощью переменной NCVLX задается число контрольных объемов в направлении координаты х, а с помощью переменной NCVLY -- соответствующее их число в направлении координаты у (для полярной системы координат переменная XL задает размер области по углу, рад). На выходе из процедуры ZGRID определены значения переменных Ll, Ml и координаты граней контрольных объемов, хранящиеся в массивах XU(I) и YV(J). Процедура ZGRID может быть использована и для введения некоторых простых неравномерностей при построении сетки. С помощью формулы

где п -- положительное число, можно получить соответствующую неравномерную сетку при различных п.

При п > 1 сетка оказывается мельче около левой границы, а вблизи правой границы она становится грубой и почти равномерной. При n < 1 сетка грубая у левой границы и мелкая и почти равномерная у правой границы.

Такое разбиение у левой и правой границ может быть изменено на противоположное, если использовать формулу

В процедуре ZGRID координаты граней XU (I) рассчитываются по (4.1) и (4.2). При п = 1, что задается по умолчанию, рассчитывается равномерная сетка. Другие значения п приводят к соответствующим неравномерностям. Значение п для оси х содержится в переменной POWERX. Хотя значения п в (4.1) и (4.2) всегда должны быть положительными, мы используем знак POWERX для выбора между формулами (4.1) и (4.2). Если POWERX > 0, EZGRID использует формулу (4.1), если же POWERX < 0, -- (4.2) (где п = ABS (POWERX)). Неравномерность разбиения по оси у создается точно так же. Для задания степени неравномерности в этом направлении используется переменная POWERY.

По умолчанию POWERX и POWERY равны единице. Если вы не измените эти значения, то получите равномерную сетку. Обычно нежелательно использовать значения POWERX и POWERY больше 2 и меньше 0,5, так как это приводит к чрезмерной неравномерности сетки.

Процедура PRINT

Когда найдено решение физической задачи, желательно вывести на печать информацию о сетке и двумерные поля зависимых переменных. Такой вывод осуществляется при вызове процедуры PRINT. Хотя обычно PRINT вызывается тогда, когда уже достигнуто сошедшееся решение, в принципе эту утилиту можно вызывать в любое время. Удобным местом для вызова PRINT является процедура OUTPUT, где желателен дополнительный вывод результатов, обеспечиваемый PRINT.

Вывод результатов процедурой PRINT включает в себя распечатку значений Х(1) и Y(J), предваряющую двумерные поля тех F(I, J,NF), для которых значение KPRINT(NF) отлично от нуля. Для каждого NF в целях идентификации поля F распечатывается заголовок TITLE (NF) длиной до 18 символов, задаваемый пользователем.

Если вы вызываете PRINT более 1 раза (например, для получения полей температуры в различные моменты времени в нестационарной задаче), то можно избежать распечатки каждый раз одних и тех же значений X (I) и Y (J). Для этой цели служит переменная KPGR. Обычно KPGR = 1, при этом вы получаете распечатку сеточных характеристик. При KPGR - 0 сетка не выводится на печать. Если вы вызываете PRINT несколько раз, то можете установить KPGR = 0 после первого вызова PRINT.

Процедура GRID

Функция процедуры GRID -- предоставление информации о расчетной сетке. В частности, необходимо задать значения MODE, LI, Ml, XU(I) для I = 2, ..., L2 и YV( J) для J = 2, ..., M2. При полярной системе координат дополнительно должно быть задано значение радиуса R(l) для нижней границы области. Грани контрольных объемов не обязательно должны равномерно отстоять одна от другой, но их координаты XU {I) и YV(J) с увеличением I и J соответственно должны расти. Для полярной системы координат значения YV

(J) должны увеличиваться в радиальном направлении. Значения XU(2) или YV (2) не обязательно должны равняться нулю. Значение R(l) может быть ненулевым, но никогда не может быть отрицательным.

Если требуется равномерная сетка, то вместо XU (I), YV (J), L1 и Ml следует задать длины XL и YL, число контрольных объемов в переменных NCVLX и NCVLY, а затем вызвать ZGRID. Можно также задать некоторые простые неравномерности с помощью значений POWERX или POWERY.

Использование неравномерной сетки является мощным средством эффективного расположения заданного числа расчетных точек. Можно легко учесть разрывы в граничных условиях, свойствах материала и распределении источников, совместив места разрывов с гранями контрольных объемов. Для создания приемлемой неравномерной сетки часто бывает полезно использование предварительных расчетов на грубых сетках, которые можно создавать с помощью ZGRID, а также самостоятельно разработав GRID.

При построении GRID помните, что ни одна составляющая неизменяемой части программы (за исключением DEFLT) еще не была выполнена. Поэтому ни геометрические величины, ни целые числа, такие как L2 или М2, еще не были рассчитаны. То, что GRID является самой первой вызываемой процедурой, предполагает множество возможностей при ее написании. Например, можно распечатать заголовок для задачи, открыть какие-нибудь файлы, выполнить некоторые предварительные действия и ввести данные оператором READ. В примерах, приведенных в книге, мы будем задавать в GRID желаемые значения текстовых переменных HEADER, PRINTF и PLOTF, которые соответствуют заголовку задачи, имени файла для вывода результатов и имени файла для графической обработки.

Процедура START

Основная функция BEGIN заключается в определении начальных значений F(I,J,NF) для соответствующих NF. Если где-либо на границе известны значения F, то желательно сразу же их задать в соответствующих граничных точках для каждого F{I,J,NF). Эти значения останутся неизменными, если соответствующие значения КВС в PHI будут сохранены равными единице.

В BEGIN удобно задавать начальные значения других переменных или параметры, остающиеся неизменными в процессе вычислений, а также большое число значений для свойств материала, плотностей тепловых потоков, граничных условий и других определяющих параметров.

Процедура OUTPUT

В то время как другие части программ должны предоставлять некоторый минимум информации для неизменяемой части программы, в отношении процедуры OUTPUT такие требования не выдвигаются. Программа будет хорошо работать, даже если оставить OUTPUT пустой. Однако (с точки зрения пользователя) OUTPUT имеет самую важную функцию. Основная задача OUTPUT -- обеспечить вывод результатов на печать. Пользователь имеет полную свободу в разработке деталей вывода результатов.

Обычно желательно вывести краткую информацию (не больше одной-двух строк) после каждой итерации и распечатать двумерные ноля нужных переменных ф после заключительной итерации. Краткий вывод после каждой итерации полезен для наблюдения за процессами сходимости стационарного решения или эволюции во времени нестационарного решения. Распечатка двумерных полей может быть произведена вызовом процедуры PRINT, созданной специально для этого.

В OUTPUT также может быть введен какой-либо критерий сходимости для остановки вычислений. Например, можно следить за изменениями некоторой представляющей интерес величины (коэффициента трения, теплового потока, максимальной температуры в области и др.) и сделать KSTOP ненулевым, когда изменения от итерации к итерации станут достаточно малыми, или же можно наблюдать за изменениями F(I, J, NF) в некоторой выбранной точке.

Так как подпрограмма OUTPUT вызывается 1 раз на каждой итерации, то в ней удобно изменять любые величины, которые меняются в зависимости от номера итерации или временных шагов. Если нужно изменить значения переменной KSOLVE(NF) для «включения» или «выключения» решения уравнения для соответствующего ф, то это должно быть сделано в OUTPUT.

Процедура РНI

В процедуре PHI содержится, наверное, наиболее важная информация о решаемой задаче. Для задач, описываемых сложными математическими моделями, PHI может быть достаточно громоздкой. В то же время для большинства часто встречающихся задач используются и довольно простые варианты PHI.

Неизменяемая часть процедуры обеспечивает решение обобщенных дифференциальных уравнений для большого числа различных переменных ф. Физический смысл F(I,J,NF) не предопределен. Как же программа узнает смысл и физическое поведение каждой зависимой переменной ф? Это происходит не с помощью информации, записанной в TITLE (NF), а с помощью данных, определенных в PHI для каждого значения NF.

Значение KSOLVE(NF) для каждого NF проверяется в подпрограмме HEART. Если KSOLVE (NF) =0, то начинается процесс расчета коэффициентов дискретных аналогов для соответствующего F (I, J, NF). В начале этого процесса происходит вызов PHI. Задача процедуры PHI заключается в определении и последующем предоставлении значений л, Г, Sc, SP и информации о граничных условиях программы. Величинам л и Г соответствуют массивы ALAM(I,J) nGAM(I,J) при 1=2, ..., L2 и J=2, ...,M2. В PHI также определяются значения коэффициентов источниковых членов Sc и Sp. Они задаются массивами SC (I, J) nSP(I,J) для 1=2, .... L2 и J=2, ..., М2. Единицы измерения Sc и Sp должны соответствовать единицам измерения скорости генерации в единице объема. Кроме того, Sp должен быть отрицательным или равным нулю. Значения Sc и Sр задаваемые по умолчанию, равны нулю, поэтому если в исходном уравнении отсутствуют источниковые члены, то не требуется выполнять никаких дополнительных действий. Наконец, задача PHI заключается в задании граничных условий с помощью индикаторов KBCIl(J), KBCLl(J), KBCJl(I), KBCMl(I). Все индикаторы КВС, заданные по умолчанию, равны единице.

5.5 Интерфейс программы

После запуска программы пользователь видит окно, в котором необходимо ввести исходные температуры, выбрать варианты вывода графической модели.

На рисунках 5.2 и 5.3 интерфейс, уже после решения задачи. Пользователю сообщается, какие параметры были применены для расчета. Можно сохранять все результаты как графические, так и табличные значения. Это поможет в дальнейшем, производить какие либо исследования уже по сохраненным результатам.

В настройках программы, предусмотрена возможность изменения конфигурации канала, вязкость жидкости, теплопроводность материала кожуха и так далее.

Также можно вывести табличные значения и для нагрева жидкости. В таблице строки соответствуют секторам, а столбцы - радиусам при построении сетки. На основе этих значений строиться графическая модель (рисунок 5.5).

6. Результаты исследования

6.1 Исследование влияния различных параметров на физические процессы

В программе, как уже упоминали раньше, предусмотрен выбор жидкости, с различной вязкостью м (Па*с); материала кожуха, у которых отличается теплопроводность k (вт/К*м); изменение конфигурации трубы, тем самым придать жесткость канала и увеличить площади теплообмена. В зависимости какие параметры были применены при расчете задачи, получаем различные температурные поля.

Приведем несколько экспериментов. Решим задачу, используя разные материалы и изменяя жидкость.

Сравним результаты вычисления для воды и для мазута. На рисунке 6.1 и 6.2 можно увидеть табличные значения нагрева этих жидкостей. В обоих случаях материал кожуха чугун.

Так же можем сравнить распределение температуры для кожуха. На рисунке 6.1 и 6.3 приведены результаты нагрева кожуха для случаев, когда материал кожуха чугун и алюминий.

6.2 Построение и объяснение графиков

Рассмотрим распределение температуры, когда конфигурация трубы разная. Выберем материал кожуха алюминий, в качестве жидкости воду.

Во всех трех случаях температура нагрева кожуха - 100 градусов, а температура жидкости изначально равно 10 градусам. По рисункам можно определить, что площадь нагрева будет больше, конечно же, у канала с шестью ребрами. Температура понижается по направлению к середине канала.

6.3 Программа расчёта и моделирования физических полей ELCUT

ELCUT® - это мощный современный комплекс программ для инженерного моделирования электромагнитных, тепловых и механических задач методом конечных элементов. Дружественный пользовательский интерфейс, простота описания даже самых сложных моделей, широкие аналитические возможности комплекса и высокая степень автоматизации всех операций позволяют разработчику полностью сосредоточиться на своей задаче.

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

ELCUT широко используется в научных исследованиях, промышленности и образовании. Академические и прикладные научные центры используют моделирование полей для изучения разнообразных физических процессов, а также для разработки приборов и устройств, применяемых в научных исследованиях.

ELCUT становится стандартным инструментом для обучения студентов различных специальностей и все более широко входит в практику преподавания физики, электротехники, прикладных дисциплин во множестве высших учебных заведений России и стран СНГ.

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

ELCUT - это программа расчёта и моделирования физических полей, и все её применения так или иначе связаны с наукой. Традиционные направления, где применяется моделирование полей, - физические исследования, электротехника, строительная механика, - в последнее время дополняются многочисленными примерами из области нанотехнологий, биологии, химии.

Моделируемые в ELCUT физические задачи, связанные с электромагнитными полями, электрическими токами, температурными полями, упругими напряжениями и деформациями - сегодня возникают в самых различных направлениях научных исследований. Большая часть применений, описанных на страницах Промышленность и Образование, тоже можно отнести к научным. Популярность программы в ведущих российских и зарубежных научных центрах доказывается её внедрениями и публикациями, связанными с ELCUT.

6.4 Проверка корректности работы программы

Результаты, полученные в процессе вычислений в дипломной работе, отвечают условиям адекватности. Чтобы убедиться в корректности работы приложения, приведем результаты, вычисленные в программе Elcut.

Температурное поле вычисляется в программе Elcut при таких же граничных условиях и температурных потоках, как и в программе реализованное нами. Если сравнить рисунок 6.7, сделанный в программном комплексе Elcut, и рисунок 6.8, построенный с помощью нашего приложения, то значения температур совпадают.

Заключение

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

Требуется исследовать, как выглядит температурное поле жидкости при изменении настраиваемых параметров.

При написании дипломной работы были использованы знания, полученные в процессе учебы.

При выполнении поставленной задачи была:

– создана математическая модель распределения температуры жидкости с заданной вязкостью в цилиндрическом канале;

– разработан алгоритм реализации модели с помощью выбранных методов;

– разработана программа по алгоритму с помощью выбранного программного обеспечения;

– разработана цифровая модель изменения поля температуры в зависимости от вязкости жидкости, теплопроводности материала кожуха, конфигурации канала и граничных условий;

– проведена проверка адекватности модели поставленной физической задачи.

В данном проекте также приведены теоретические сведения, необходимые для решения поставленной задачи и рассмотрены результаты выполненной работы.

Все, поставленные ранее задачи, выполнены. Полученный программный продукт готов к использованию.

Список использованной литературы

1. Патанкар С.В., Численное решение задач теплопроводности и конвективного теплообмена при течениях в каналах: пер. с англ. Е.В. Калабина; под ред Г.Г. Янькова. - М.: Издательство МЭИ, 2003. - 312с., ил

2. Чижиумов С.Д., Основы гидродинамики: учеб. пособие - Комсомольск-на-Амуре : ГОУВПО «КнАГТУ», 2007. - 106 с.

3. Патанкар С.В., Численные методы решения задач теплообмена и динамики жидкости: пер. с англ. - М. : Энергоатомиздат, 1984. - 152 с.

4. Сологаев В.И. Гидравлика (механика жидкости и газа): учеб. пособие. - Омск: изд. СибАДИ, 2010. - 64с.

5. Выгодский М.Я. Справочник по элементарной математике.- М.: Наука, 1982.- 336с.

6. Бронштейн И.Н., Семендяев К.А. Справочник по математике для инженеров и учащихся ВТУЗов.- 13-е изд., исправленное.-М.: Наука, 1986.- 544с.

7. Елизарова Т.Г., Математические модели и численные методы в динамике газа и жидкости. - М.: Физический факультет МГУ, 2005. - 224с.

8. Гебхарт Б., Джалурия Й., Свободноконветивные течения, тепло- и массообмен., кн. 1., пер. с англ. - М.: Мир, 1991. - 678 с., ил.

Приложение

procedure DEFLT // задание переменных по умолчанию

var

nz,n,i,j:integer;

begin

kstop:=0;

last:=5;

iter:=0;

kord:=2;

mode:=1;

kpgr:=1;

kout:=3;

small:=1.e-20;

big:=1.e+20;

time:=0.;

dt:=1.e+20;

r[1]:=0.;

powerX:=1.;

powerY:=1.;

for nz:=1 to nzmax do

begin

powrX[nz]:=1.;

powrY[nz]:=1.;

end;

end;

procedure PHI // определение значений

procedure Phi;

var

i,j:integer;

begin

for j:=2 to M2 do

for i:=2 to L2 do

begin

gam[i,j]:=cond; //коэфф-т диффузии Г

end;

for j:=2 to M2 do // Определение граничных условий

begin

kbci1[j]:=2;

kbcl1[j]:=2;

end;

for i:=1 to ncvlx do // Определение граничных условий

// для представления правой границы в виде оси симметрии

kbcm1[i]:=2;

end;

procedure GRID // предоставление информации о расчетной сетке

procedure Grid;

begin

rout:=1.; // радиус кожуха

rin:=0.8; // радиус трубы

alpha:=180; // рассчитываемая область, половина круга

xl:=0.5*alpha*Pi*rout*rout/180; // формула круга

r[1]:=rin; // радиус вектор

yl:=rout-rin; // толщина кожуха

ncvlx:=47; // количество рассчитываемых радиус векторов

ncvly:=15; // количество рассчитываемых секторов

mode:=3; // применяется полярная система координат

EZgrid; // вызов процедуры

end;

procedure BEGIN // Определение начальных граничных условий

var

i,j:integer;

begin

ksolve[2]:=1;

kprint[1]:=1;

kprint[2]:=1;

kplot[1]:=1;

kplot[2]:=1;

last:=6;

// задаем теплопроводность металлам

if Form1.N14.Checked then cond:=3;

if Form1.N15.Checked then cond:=7;

if Form1.N16.Checked then cond:=11;

if Form1.N17.Checked then cond:=19;

//присваиваем произвольные значения

cp:=1; //теплоемкости

den:=1; //плотности

dpdz:=-1; //градиенту давления

t1:=StrToFloat(Form1.Edit1.Text); // температура нагрева

t2:=StrToFloat(Form1.Edit2.Text); // температура стенки канала

temperGraf := t2;

phocp:=den*cp; //объемная теплоемкость

for j:=1 to m1 do

for i:=1 to l1 do

begin

f[i,j,1]:=0;

f[i,j,2]:= t1; //t1;

f[i,m1,2]:= t2; //t2;

end;

end;

procedure OUTPUT; // расчет и вывод на экран

var

j,i,iunit:integer;

begin

if(iter=0) then

begin

ksolve[1]:=1;

end;

if(iter=3) then // после трех итераций нахождения скорости w,

begin // начинаем находить температуру Т

ksolve[1]:=1;

ksolve[2]:=0;

end;

Print;

end;

procedure PRINT // Вывод информации о сетке и двумерных полей

var

i,j:integer;

temper, temper1: double;

begin

// выводим значения в таблицу

for i:=49 downto 1 do

for j:=1 to 17 do

Form3.StringGrid2.Cells[i,j] := FloaToStrF((f[i,j,2]),ffFixed,5,2);

// выводим заголовки таблицы

for i:=49 downto 1 do

Form3.StringGrid2.Cells[i,0] := IntToStr(i) + ' радиус';

for j:=1 to 17 do

Form3.StringGrid2.Cells[0,j] := IntToStr(j) + ' сектор';

temper := (StrToFloat(Form3.StringGrid2.Cells[48,16]) +

StrTFloat(Form3.StringGrid2.Cells[49,16]))/2;

Form3.StringGrid2.Cells[48,17]:= FloatToStrF(temper,ffFixed,5,2);

temper := (StrToFloat(Form3.StringGrid2.Cells[2,16]) +

StrTFloat(Form3.StringGrid2.Cells[1,16]))/2;

Form3.StringGrid2.Cells[1,17]:= FloatToStrF(temper,ffFixed,5,2);

trnagrev := StrToFloat(Form3.StringGrid2.Cells[13,2]);

trnagrev2 := StrToFloat(Form3.StringGrid2.Cells[1,2]);

if (Form1.N41.Checked = true) then

trnagrev1 := StrToFloat(Form3.StringGrid2.Cells[7,2])

else

if (Form1.N61.Checked = true) then

begin

trnagrev1 := StrToFloat(Form3.StringGrid2.Cells[4,2]);

trnagrev3 := StrToFloat(Form3.StringGrid2.Cells[9,2]);

end;

end;

procedure ZGRID // предоставление информации о расчетной сетке

var

i,i1,i2,j,j1,j2,nz,ilast,jlast:integer;

dd,fcvlx,fcvly:double;

begin

xu[2]:=0.;

i2:=2;

for nz:=1 to nzx do

begin

fcvlx:=(ncvx[nz]);

ilast:=i2;

i1:=ilast+1;

i2:=ilast+ncvx[nz];

for i:=i1 to i2 do

begin

dd:=(i-ilast)/fcvlx;

if (powrx[nz]>0) then

xu[i]:=xu[ilast]+xzone[nz]*power(dd,powrx[nz])

else

xu[i]:=xu[ilast]+xzone[nz]*(1.-power((1.-dd),(-powrx[nz])));

end;

end;

L1:=i2; //Рассматривается направление Y

yv[2]:=0.;

j2:=2;

for nz:=1 to nzy do

begin

fcvly:=ncvy[nz];

jlast :=j2;

j1:=jlast+1;

j2:=jlast+ncvy[nz];

for j:=j1 to j2 do

begin

dd:=(j-jlast)/fcvly;

if (powry[nz]>0) then

yv[j]:=yv[jlast]+yzone[nz]*power(dd,powry[nz])

else

yv[j]:=yv[jlast]+yzone[nz]*(1.-power((1.-dd),(-powry[nz])));

end;

end;

M1:=j2;

end;

procedure EZGRID

var

fcvlx,fcvly,dd:double;

i,j:integer;

begin

//Построение сетки зона за зоной

//Рассматривается направление X

L1:=ncvlx+2;

xu[2]:=0.; //координаты граней по х

xu[L1]:=xl;

L2:=L1-1;

fcvlx:=ncvlx;

for i:=3 to L2 do

begin

dd:=(i-2)/fcvlx;

if (powerx>0) then

xu[i]:=xl*power(dd,powerx)

else

xu[i]:=xl*(1.-power((1.-dd),(-powerx)));

end;

M1:=ncvly+2;

yv[2]:=0.;

yv[M1]:=yl;

M2:=M1-1;

fcvly:=ncvly;

for j:=3 to M2 do

begin

dd:=(j-2)/fcvly;

if (powery>0) then

yv[j]:=yl*power(dd,powery)

else

yv[j]:=yl*(1.-power((1.-dd),(-powery)));

end;

end;

procedure HEART()

// Расчет коэффициентов дискретных аналогов и модификация граничных условий

label

Lab2;

var

n,j,i:integer;

beta,rlx,vol,apt,diff,temp,area,anb,ainr:double;

begin

//Построение цикла для всех аналогов

for n:=1 to nfmax do //последовательное вычисление каждой зависимой переменной

begin

nf:=n;

if (ksolve[nf]=0) then goto Lab2; //переход к следующей итерации цикла

Phi;

//Вычисление коэффициентов дискретных аналогов

beta:=4./3.;

if (kord=1) then beta:=1.;

rlx:=(1.-relax[nf])/relax[nf];

//Рассматриваются объемные условия

for j:=2 to M2 do

for i:=2 to L2 do

begin

vol:=Ycvr[j]*Xcv[i]; //объем КО

apt:=alam[i,j]/dt;

sc[i,j]:=(sc[i,j]+apt*f[i,j,nf])*vol;//коэф-т b

sp[i,j]:=(apt-sp[i,j])*vol; //коэфф-т ар

end;

//коэффициенты для распространения по Х

for j:=2 to M2 do

for i:=2 to L3 do

begin

diff:=arx[j]*2.*gam[i,j]*gam[i+1,j]/((Xcv[i]*gam[i+1,j]+Xcv[i+1]*gam[i,j]+small)*sx[j]);

aip[i,j]:=diff+small; //коэф-т aE

aim[i+1,j]:=aip[i,j]; //коэф-т aW

end;

for j:=2 to M2 do

begin

//Рассматривается граница i=1

diff:=gam[2,j]/(0.5*xcv[2]*sx[j])+small; //диффузионная

aim[2,j]:=beta*diff; //коэф-т aW

aip[1,j]:=aim[2,j]; //коэф-т aE

aim[2,j]:=aim[2,j]*arx[j];

aim[1,j]:=(beta-1)*aip[2,j]/arx[j];

aip[2,j]:=aip[2,j]+aim[1,j]*arx[j];

if (kbci1[j]=1) then

sc[2,j]:=sc[2,j]+aim[2,j]*f[1,j,nf]

else

begin

sp[1,j]:=aip[1,j]-flxpi1[j]; //коэф-т aР

sc[1,j]:=flxci1[j]; //коэф-т b

temp:=aim[2,j]/sp[1,j];

sp[2,j]:=sp[2,j]-aip[1,j]*temp;

aip[2,j]:=aip[2,j]-aim[1,j]*temp; //коэф-т aE

sc[2,j]:=sc[2,j]+sc[1,j]*temp;

end;

sp[2,j]:=sp[2,j]+aim[2,j];

aim[2,j]:=0.;

//Рассматривается граница i=L1

diff:=gam[L2,j]/(0.5*xcv[L2]*sx[j])+small;

aip[L2,j]:=beta*diff;

aim[L1,j]:=aip[L2,j];

aip[L2,j]:=aip[L2,j]*arx[j];

aip[L1,j]:=(beta-1)*aim[L2,j]/arx[j];

aim[L2,j]:=aim[L2,j]+aip[L1,j]*arx[j];

if (kbcL1[j]=1) then

sc[L2,j]:=sc[L2,j]+aip[L2,j]*f[L1,j,1]

else

begin

sp[L1,j]:=aim[L1,j]-flxpL1[j];

sc[L1,j]:=flxcL1[j];

temp:=aip[L2,j]/sp[L1,j];

sp[L2,j]:=sp[L2,j]-aim[L1,j]*temp;

aim[L2,j]:=aim[L2,j]-aip[L1,j]*temp;

sc[L2,j]:=sc[L2,j]+sc[L1,j]*temp;

end;

sp[L2,j]:=sp[L2,j]+aip[L2,j];

aip[L2,j]:=0.;

end;

//Коэффициенты для распространения Y

for j:=2 to M3 do

for i:=2 to L2 do

begin

area:=rv[j+1]*xcv[i]; //площадь грани КО

diff:=area*2.*gam[i,j]*gam[i,j+1]/(ycv[j]*gam[i,j+1]+ycv[j+1]*gam[i,j]+small);

ajp[i,j]:=diff+small; //коэф-т aN

ajm[i,j+1]:=ajp[i,j]; //коэф-т aS

end;

for i:=2 to L2 do

begin

//Рассматривается граница j=1

area:=rv[2]*xcv[i];

diff:=gam[i,2]/(0.5*ycv[2])+small;

ajm[i,2]:=beta*diff;

ajp[i,1]:=ajm[i,2];

ajm[i,2]:=ajm[i,2]*area;

ajm[i,1]:=(beta-1)*ajp[i,2]/(rv[3]*xcv[i]);

ajp[i,2]:=ajp[i,2]+ajm[i,1]*area;

if (kbcj1[i]=1) then

sc[i,2]:=sc[i,2]+ajm[i,2]*f[i,1,nf]

else

begin

sp[i,1]:=ajp[i,1]-flxpj1[i];

sc[i,1]:=flxcj1[i];

temp:=ajm[i,2]/sp[i,1];

sp[i,2]:=sp[i,2]-ajp[i,1]*temp;

ajp[i,2]:=ajp[i,2]-ajm[i,1]*temp;

sc[i,2]:=sc[i,2]+sc[i,1]*temp;

end;

sp[i,2]:=sp[i,2]+ajm[i,2];

ajm[i,2]:=0.;

//Рассматривается граница j=M1

area:=rv[M1]*xcv[i];

diff:=gam[i,M2]/(0.5*ycv[M2])+small;

ajp[i,M2]:=beta*diff;

ajm[i,M1]:=ajp[i,M2];

ajp[i,M2]:=ajp[i,M2]*area;

ajp[i,M1]:=(beta-1)*ajm[i,M2]/(rv[M2]*xcv[i]);

ajm[i,M2]:=ajm[i,M2]+ajp[i,M1]*area;

if (kbcM1[i]=1) then

sc[i,M2]:=sc[i,M2]+ajp[i,M2]*f[i,M1,2]

else

begin

sp[i,M1]:=ajm[i,M1]-flxpM1[i];

sc[i,M1]:=flxcM1[i];

temp:=ajp[i,M2]/sp[i,M1];

sp[i,M2]:=sp[i,M2]-ajm[i,M1]*temp;

ajm[i,M2]:=ajm[i,M2]-ajp[i,M1]*temp;

sc[i,M2]:=sc[i,M2]+sc[i,M1]*temp;

end;

sp[i,M2]:=sp[i,M2]+ajp[i,M2];

ajp[i,M2]:=0.;

end;

//Построение массивов sp[i,j] и sc[i,j]

for j:=2 to M2 do

begin

for i:=2 to L2 do

begin

anb:=aip[i,j]+aim[i,j]+ajp[i,j]+ajm[i,j]; //сумма соседних коэф-тов

ainr:=anb*rlx; //инерция I

sp[i,j]:=sp[i,j]+anb+ainr; //коэф-т aP

sc[i,j]:=sc[i,j]+ainr*f[i,j,1]; //коэф-т b

end;

//Вызов процедуры Solve для получения решения дискретных аналогов

end;

Solve;

end;

Lab2:

time:=time+dt;

iter:=iter+1;

if (iter>=last) then kstop:=1;

end;

procedure TForm1.Button1Click(Sender: TObject);

label

Lab1;

begin

if not TryStrToFloat(Form1.Edit1.Text, t1) or not TryStrT

Float(Form1.Edit2.Text, t2)

or (StrToFloat(Form1.Edit2.Text) < 0) or (StrTFloat(Form1.Edit1.Text) < 0) then

begin

MessageBox(handle, PChar('Введите корректные температуры!'), PChar('Ошибка'), MB_OK+MB_ICONWARNING)

end

else

begin

if t1 >= t2 then

begin

MessageBox(handle, PChar(' Температура на границах внутренней трубы должна ' + #10 + 'быть не равной и меньше температуры на участке нагрева!'), PChar('Ошибка'), MB_OK+MB_ICONSTOP)

end

else

begin

if StrToFloat(Form1.Edit3.Text) >= StrToFloat(Form1.Edit1.Text) then

begin

MessageBox(handle, PChar('Температура жидкости должна быть меньше температуры стенки канала!'), PChar('Ошибка'), MB_OK+MB_ICONSTOP)

end

else

if (StrToFloat(Form1.Edit2.Text) - StrToFloat(Form1.Edit1.Text) < 10) then

begin

end

else

begin

Form1.Label4.Caption := 'Выбран материал кожуха: ' + material + #10 + 'На внутренней границе ' + Form1.Edit1.Text + ' C' + #176 + #10 + 'По каналу течет жидкость ' + Form1.Edit3.Text + ' C' + #176 + #10 + 'Нагрев идет с нижней части трубы ' + Form1.Edit2.Text + ' C' + #176;

if (Form1.RadioGroup1.ItemIndex = 2) then

begin

Form3.Label1.Caption := 'Материал кожуха: ' + material + '.';

Form3.Caption := 'Нагрев кожуха ';

end;

if (Form1.RadioGroup1.ItemIndex = 1) then

begin

Form3.Label1.Caption := 'Выбрана жидкость:' + gidkost + '.';

Form3.Caption := 'Нагрев жидкости';

end;

Deflt;

Grid;

Ready;

Start;

Lab1:

Output;

if (kstop<>0) then

begin

metka := 0;

Form1.Button3.Click;

abort;

end;

Heart;

goto Lab1;

abort;

end;

end;

end;

end;

procedure TForm1.Button2Click(Sender: TObject);

// нажатие на кнопку строит графическую иллюстрацию

var

i, j: integer;

begin

if Form3.StringGrid2.Cells[1,1] = '' then

MessageBox(handle, PChar('Сначала необходимо произвести расчет! !'), PChar('Ошибка'), MB_OK+MB_ICONERROR)

else

begin

Form4.N1.Click;

Form4.Show();

end;

end;

procedure TForm3.Button1Click(Sender: TObject);

// сохранение значений в excel файл

var

znac : double;

begin

if (Form1.RadioGroup1.ItemIndex = 2) then

begin

if (Form3.stringGrid2.Cells[49,17] = '') then

MessageBox(handle, PChar('Таблица пустая! ' + #10 + 'Вы точно хотите сохранить?'), PChar('Сохранение'), MB_OK+MB_ICONWARNING)

else

if SaveAsExcelFile(Form3.stringGrid2, 'Значения нагрева кожуха', 'c:teplo.xls') then

ShowMessage('Данные успешно сохранены!');

end;

if (Form1.RadioGroup1.ItemIndex = 1) then

begin

if (Form3.stringGrid1.Cells[1,17] = '') then

MessageBox(handle, PChar('Таблица пустая! ' + #10 + 'Вы точно хотите сохранить?'), PChar('Сохранение'), MB_OK+MB_ICONWARNING)

else

if SaveAsExcelFile(Form3.stringGrid1, 'Значения нагрева жидкости', 'c:teplo1.xls') then

ShowMessage('Данные успешно сохранены!');

end

end;

procedure TForm4.Button1Click(Sender: TObject);

// построение графической иллюстрации

var

i,j1,l,l1,mass,mass1 :integer;

temp, tem1,znach : double;

begin

// задается стиль пера

Form4.image1.Canvas.Brush.Style := bsDiagCross;

// задается цвет карандаша

Form4.image1.Canvas.Pen.Color := clBlack;

// строиться прямоугольник

Form4.image1.Canvas.Rectangle(0,0,800,650);

// очищается стиль

Form4.image1.Canvas.Brush.Style := bsClear;

tem1 := (temperGraf+2)/210;

l := 38;

l1 := l;

j := 7;

// строиться сам рисунок, для каждой четверти круга свой цикл

// левый верхний

for mass:=26 to ncvlx+2 do

begin

j1 := j;

mass1 :=75 - mass;

for i:= ncvly+2 downto 1 do

begin

// сравнивается код цвета со значением и красится в определеннй

// цвет

Form4.image1.Canvas.Brush.Color := color_schema[round(StrToFloat(Form3.StringGrid2.Cells[mass1,i])/tem1)];

Form4.image1.Canvas.Pen.Color := clBlack;

if Form4.N2.Checked = false then

Form4.image1.Canvas.Pen.Color := color_schema[round(StrToFloat(Form3.StringGrid2.Cells[mass1,i])/tem1)];

Form4.image1.Canvas.Pie(j1,j1,650-j1,650-j1,massI[mass+1],massJ[mass+1],massI[mass],massJ[mass]);

j1 := j1+j;

end;

end;

// левый нижний

for mass:=1 to 24 do

begin

j1 := j;

mass1 := mass + 1;

for i:=17 downto 1 do

begin

Form4.image1.Canvas.Brush.Color := color_schema[round(StrToFloat(Form3.StringGrid2.Cells[25-mass,i])/tem1)];

Form4.image1.Canvas.Pen.Color := clBlack;

if Form4.N2.Checked = false then

Form4.image1.Canvas.Pen.Color := color_schema[round(StrToFloat(Form3.StringGrid2.Cells[25-mass,i])/tem1)];

Form4.image1.Canvas.Pie(j1,j1,650-j1,650-j1,massI[mass+1],massJ[mass+1],massI[mass],massJ[mass]);

j1 := j1+j;

end;

end;

// правый верхний

for mass:=1 to 24 do

begin

j1 := j;

mass1 := 25 - mass;

for i:=17 downto 1 do

begin

Form4.image1.Canvas.Brush.Color := color_schema[round(StrToFloat(Form3.StringGrid2.Cells[mass1,i])/tem1)];

Form4.image1.Canvas.Pen.Color := clBlack;

if Form4.N2.Checked = false then

Form4.image1.Canvas.Pen.Color := color_schema[round(StrToFloat(Form3.StringGrid2.Cells[mass1,i])/tem1)];

Form4.image1.Canvas.Pie(j1,j1,650-j1,650-j1,massI1[mass],massJ1[mass],massI1[mass+1],massJ1[mass+1]);

j1 := j1+j;

end;

end;

// правый нижний

for mass:=26 to 49 do

begin

j1 := j;

mass1 := 75 - mass;

for i:=17 downto 1 do

begin

Form4.image1.Canvas.Brush.Color := color_schema[round(StrToFloat(Form3.StringGrid2.Cells[mass1,i])/tem1)];

Form4.image1.Canvas.Pen.Color := clBlack;

if Form4.N2.Checked = false then

Form4.image1.Canvas.Pen.Color := color_schema[round(StrToFloat(Form3.StringGrid2.Cells[mass1,i])/tem1)];

Form4.image1.Canvas.Pie(j1,j1,650-j1,650-j1,massI1[mass],massJ1[mass],massI1[mass+1],massJ1[mass+1]);

j1 := j1+j;

end;

end;

Form4.image1.Canvas.Brush.Style := bsClear;

Form4.image1.Canvas.Brush.Style := bsSolid;

Form4.image1.Canvas.Brush.Color := clWhite;

Form4.image1.Canvas.Pen.Color := clBlack;

Form4.image1.Canvas.Ellipse(130,130,520,520);

Form4.image1.Canvas.Brush.Style := bsClear;

end;

ref.by 2006—2019
contextus@mail.ru