Научная статья на тему 'Дискретизация дифференциального уравнения конвекции и диффузии на основе метода контрольного объема'

Дискретизация дифференциального уравнения конвекции и диффузии на основе метода контрольного объема Текст научной статьи по специальности «Математика»

CC BY
457
59
i Надоели баннеры? Вы всегда можете отключить рекламу.
Ключевые слова
CONVECTION / HEAT AND MASS TRANSFER / MATHEMATICAL MODEL / GENERALIZED DIFFERENTIAL EQUATION / DISCRETE ANALOGUE / CONTROL VOLUME / FUNCTION APPROXIMATION / КОНВЕКЦИЯ / ТЕПЛОМАССОПЕРЕНОС / МАТЕМАТИЧЕСКАЯ МОДЕЛЬ / ОБОБЩЕННОЕ ДИФФЕРЕНЦИАЛЬНОЕ УРАВНЕНИЕ / ДИСКРЕТНЫЙ АНАЛОГ / КОНТРОЛЬНЫЙ ОБЪЕМ / АППРОКСИМАЦИЯ ФУНКЦИИ

Аннотация научной статьи по математике, автор научной работы — Боков Александр Викторович, Клячин Алексей Александрович, Корытова Марина Александровна

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

i Надоели баннеры? Вы всегда можете отключить рекламу.

Похожие темы научных работ по математике , автор научной работы — Боков Александр Викторович, Клячин Алексей Александрович, Корытова Марина Александровна

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

DISCRETIZATION OF DIFFERENTIAL EQUATIONS OF CONVECTION AND DIFFUSION BASED ON CONTROL VOLUME METHOD

The object of research is a mathematical model of convective heat and mass transfer. This model is based on the equations of hydrodynamics. Taking into account that the equations governing heat transfer, mass transfer and fluid flow are similar in many respects, and the dependent variables of our interest satisfy a generalized conservation law, we apply numerical methods to a differential equation of a generalized variable. The research objective is construction of a discrete analogue of the generalized differential equation describing the proposed mathematical model of convection in a viscous incompressible fluid. To implement the mathematical model, it is suggested to use the control volume method. This method has a number of significant advantages in comparison to the method of finite differences. It provides strict compliance with the laws of mass, momentum, energy conservation for any group of cells, and therefore, in the whole computational domain. The method itself and its application are sufficiently described in the book by S. Patankar [5]. The computational domain is divided into a multiplicity of control volumes. The differential equation is integrated over the control volumes. In this case, the assumption about the form of the function describing the change in the variable between two adjacent nodes is made. To represent the profiles of functions between the nodes of the grid, we use templates that approximate the exact solution obtained analytically. As a rule, researchers use discrete analogues obtained for rectangular grids in numerical simulation. Use of such templates for curvilinear coordinate systems leads to large computational errors. Many applied problems are easier to solve in the curvilinear coordinates, for example, in cylindrical coordinates. Further, to obtain the exact templates, we need different solutions that take into account specific features of curvilinear grids. Modern computers allow performance of complex mathematical calculations. It also makes possible the use of more accurate functional dependencies for approximations. This all led to the desire of authors to get a more efficient tool the calculated schemes and algorithms that use an exponential scheme for the cylindrical coordinate system as a basis. The task is to obtain the best approximations of the generalized variable profiles in a variety of directions. For this purpose, we find the exact solutions of the equations of conservation separately for each coordinate. The article investigates the meaning of exact solutions made in order to ensure their correctness. Using the obtained analytical solutions, we have built the discrete analogue of the generalized differential equation. We obtained a more general formulation for the discrete analogue based on the research of S. Patankar. The resulting discrete analogue of the generalized differential equation written in cylindrical coordinates can be used for numerical simulations of convective heat and mass transfer under conditions of axisymmetric laminar and turbulent flows of viscous incompressible fluid and fluid flows in channels and pipes. Using discrete analogues constructed on the basis of exact transfer equations may be more advantageous taking into account the present level of computer technology.

Текст научной работы на тему «Дискретизация дифференциального уравнения конвекции и диффузии на основе метода контрольного объема»

www.volsu.ru

DOI: https://doi.oгg/10.15688/j•volsu1.2016.4.2

УДК 519.63 ББК 22.311

ДИСКРЕТИЗАЦИЯ ДИФФЕРЕНЦИАЛЬНОГО УРАВНЕНИЯ КОНВЕКЦИИ И ДИФФУЗИИ НА ОСНОВЕ МЕТОДА КОНТРОЛЬНОГО ОБЪЕМА

Александр Викторович Боков

Кандидат физико-математических наук,

доцент кафедры естественнонаучных и гуманитарных дисциплин, Российский экономический университет им. Г.В. Плеханова, филиал в г. Пятигорске Ставропольского края av_bokov@mail.ru

ул. Кучуры, 8, 357500 г. Пятигорск, Российская Федерация

Алексей Александрович Клячин

Доктор физико-математических наук,

заведующий кафедрой математического анализа и теории функций, Волгоградский государственный университет klyachin-aa@yandex.ru, matf@volsu.ru

просп. Университетский, 100, 400062 г. Волгоград, Российская Федерация

ю Марина Александровна Корытова

20 Кандидат физико-математических наук,

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

Южно-Уральский государственный университет

и (национальный исследовательский университет)

° sma.59@mail.ru

& просп. Ленина, 76, 454080 г. Челябинск, Российская Федерация

¡2

^ Аннотация. Объектом исследования является математическая модель

® конвективного тепломассопереноса. Цель исследования — построение дис-

¡г

« кретного аналога для обобщенного дифференциального уравнения, описыва-

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

оз мои жидкости.

т Ключевые слова: конвекция, тепломассоперенос, математическая мо-

о дель, обобщенное дифференциальное уравнение, дискретный аналог, кон-

^ трольный объем, аппроксимация функции.

Введение

Авторы статьи некоторое время занимались численным моделированием конвективного переноса. Для реализации математической модели было предложено использовать метод контрольного объема [5; 6]. Этот метод имеет ряд существенных преимуществ по сравнению с методом конечных разностей. Он гарантирует точное выполнение законов сохранения массы, количества движения, энергии на любой группе ячеек — контрольных объемов, а потому и во всей расчетной области. Для представления профилей функций между узлами сетки использовались различные шаблоны, аппроксимирующие точное решение, полученное аналитически для элементарного отрезка между соседними узлами. Применение шаблона на основе точного решения приводило к экспоненциальной схеме, требующей больших вычислительных ресурсов. Поэтому использование ее ранее считалось неэффективным. Кроме того, решения были получены для прямоугольных сеток. Дискретные аналоги и примеры расчетов на таких сетках можно найти в работах [1; 4; 7]. Использование таких же шаблонов для криволинейных систем координат приводило к большим вычислительным погрешностям, что с математической точки зрения было совершенно не оправдано. Многие прикладные задачи принципиально удобнее решать в криволинейных, например, в цилиндрических координатах. Но тогда для получения точных шаблонов нужны другие точные решения, учитывающие особенности криволинейных сеток. Современные требования к точности численных решений также изменились, и вычислительные мощности современных ЭВМ позволяют эффективно проводить сложные вычисления. Все это объясняет желание авторов получить более эффективный инструмент — расчетные схемы и алгоритмы, использующие в своей основе именно экспоненциальную схему для цилиндрической системы координат.

1. Математическая модель конвективного теплопереноса и обобщенное дифференциальное уравнение

Математическая модель свободной или вынужденной конвекции в вязкой несжимаемой жидкости получается на основе уравнений гидродинамики, тепло- и массообмена (см., например, [2]): уравнения Навье — Стокса для несжимаемой жидкости

/дй \ 1

р( д^ + (й ■ У)й) = -Ур + + (У-й) + рд, (1)

уравнения неразрывности

др + V-(рй) = 0 (2)

и уравнения переноса тепла (аналогичный вид имеет уравнение переноса массовой концентрации химической компоненты)

ср(|T + и •VT^j = AV2T. (3)

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

Дифференциальные уравнения (1)-(3) описывают процессы переноса количества движения, массы и энергии. Каждое из уравнений выражает соответствующий закон сохранения (импульса, массы, энергии). Рассматривая совместно (1-2), получаем систему уравнений для компонент вектора скорости в виде

д д'р

— (рщ) + V ■ (рииг) = -+ V ■ (цЧщ) + Вг + Ц, (4)

где щ — г-я компонента вектора скорости (г = 1, 2, 3); В^ — г-я составляющая объемной силы (приложенной к единице объема); Ц — г-я составляющая величины, учитывающей дополнительные диссипативные члены.

Из (3) и (4) следует, что щ (г = 1,2,3) и Т подчиняются обобщенному закону сохранения, которому соответствует дифференциальное уравнение для зависимой (обобщенной) переменной Ф :

д

— (рФ) + V- (риФ) = V- (^Ф) + Ф, (5)

где Г — коэффициент «диффузии», а Ф — так называемый «источниковьш член» (обычно представляется в линеаризованном виде, как линейная функция Ф). Здесь физический смысл Ф и Г зависит от того, какую именно величину обозначает переменная Ф. Коэффициент Г играет роль коэффициента диффузии для уравнения переноса массовой концентрации химической компоненты, коэффициента теплопроводности для уравнения тепловой конвекции, коэффициента динамической или кинематической вязкости для уравнения движения жидкости. Источниковый член Ф может включать диссипативные и объемные компоненты, которые не учитываются конвективным и диффузионным членами уравнения (5).

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

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

В [5] приведены примеры дискретных аналогов для уравнения (5), записанного в декартовых координатах. Их различие обусловлено способом выбора интерполяционных

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

2. Интегральные балансы для граней контрольного объема

Уравнение (5) с использованием уравнения неразрывности (2) можно упростить, представив в виде:

д Ф

р — + ри • УФ = V • (ГУФ) + Ф,

Запишем последнее уравнение в цилиндрических координатах (г, ф, г) (с учетом правил действий с оператором Гамильтона, [3]):

д Ф д Ф 1 дФ д Ф

>— + + р-т;— uv + р—1

dt дг г д ф dz

1 д ( ^ дФ\ 1 д (Г дФ\ 1 д ( ^ дФ\ т (6)

=--г Г— +-----+---г Г— + Ф.

д д дф дф д

Здесь иг,иф,их представляют собой соответствующие компоненты вектора скорости: и = [иг,иф,иг).

Преобразуем уравнение (6), умножив обе части на г :

Ф Ф дФ дФ

гр— + григ— + риф— + гри^ — дt дг дф дг

д дФ Г дФ д Ф

— г Г— +----+--г Г— + г Ф.

д д дф дф д

В цилиндрических координатах уравнение неразрывности имеет вид:

(7)

др 1 д , , 1 д , , д , , дь + гд~г ^ + гдф(ри^ + Тг [Ри) = °

или

д д д д М + д~г + + Ж= °. (8)

Преобразуем уравнение конвекции и диффузии (7), используя уравнение неразрывности (8):

д д д д

Ш( ГРФ) + Гг(гри,.ф) + дф(р..фф) + д-Ми.Ф) =

= д (гГ?*) + ± (Г »и д (ГГ?*)+г*. (9)

д д дф дф д д

Уравнение (9) проинтегрируем по контрольному объему V (рис. 1) и по временному промежутку [¿о, ¿о + Д£]. Обозначим через и ид, верхнюю и нижнюю грани контрольного

объема, через w и е — левую («западную») и правую («восточную») грани, через п и ^ — «северную» и «южную» грани.

и

Рис. 1. Контрольный объем V с узловой точкой Р

Для удобства представления формул используем такие же обозначения в записи соответствующих этим граням значений переменных. Таким образом, через обозначим значение переменной г на нижней грани, а через ги — значение этой переменной на верхней грани. Тогда для контрольного объема V выполнено условие < ^ < ги. Аналогичным образом вводятся обозначения для значений переменных ф иг, тогда фе < ф < фт, Гз < т < гп. Следовательно, интегрирование по рассматриваемому контрольному объему — это интегрирование по области

V = {(Г, ф, г) | Г3 <Г< Гп, фе < ф < фт, <х< г,и}.

Интегрируя по области V, получаем

¿О+Д*

ЯГ

V

- (грФ)^

¿г ¿ф ¿X +

¿О+Д* фw % и

- ^ (гГ I I ¿г

ог \ иг

) I ¿—хо ууы^и ' п / О

I Я ДI ф

ь + Iw / ^ ^ /'Г д Ф \ \

+ / Я Д^(р""ф) - ^Ьд*))^

е

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

К Тг <гри-Ф -ш{гГ Ш'Ь

¿0 ГаХЛ фе ¿О+Д* Т п фw

+ I я

¿0 ге фе

¿г ¿ф сИ +

¿г ¿г ей + ¿ ¿ ¿ =

¿о+Д* г

/ лгг ^ ¿г

¡о I. V

а.

(10)

Размеры контрольного объема V определяются значениями Дг = гп — г3, Дф = = фе — фт (рис. 2) и Дг = хи — га- Пусть (Ьг)п — расстояние между узловыми точками Р и N, а (6 г)3 — расстояние между Р и Б .По такому же принципу обозначим расстояния между другими парами узлов: (6ф)ад — между Р и Ш, (6ф)е — между Р и Е, (Ьг)и — между Р и и, (6г)а — между Р и И.

Ar

(Sri

n

\ \ ______ / P 1

\\ r \ \ _ / V 1 1 Аф "7 +-—-JSH

E

S

Рис. 2. Проекция контрольного объема V на плоскость (г, ф)

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

д Ф

Г дФ

дФ

Jr = rpurФ - гГ—, Jv = рМфФ--—, Jz = rpuzФ - гГ —

дг ' г д ф'

С учетом формул (11) уравнение (10) примет вид:

д

(11)

io+Ai

+

ЯГ

V

Т' ггг (_JL + __J^ + дЬ

J JJJ V _г дф dz

io+Ai

Jd

- {г РФ) dt ^ dr dф dz

dr dф dz +

io+Ai

t o

dt = J JJJ r^dr dф dz

L У

(12)

dt.

Пусть ,1п — значение потока ,1Г на грани п, а — значение ,1Г на грани 8. Также обозначим через Jw и значения ,]ф на гранях т и е, а через Зи и — значения на гранях и и d. Значения переменных Ф и р в момент времени to + Дt (новые значения на текущем временном шаге) в узловой точке Р будем обозначать как Фр и рР. Подобным образом обозначим значения переменной Ф в соседних узлах: Ф^ в узле N, Ф^ в узле Б, Ф^ и Фр в точках Ш и Е, Ф^ и Фр в точках и и И. Значения Ф и р в узле Р в момент времени ¿0 («старые», то есть заданные значения) обозначим через Фр и рр. Согласно применяемому методу считаем, что значения переменных Ф и р в узловой точке сохраняются во всем окружающем эту точку контрольном объеме. Кроме того, при интегрировании будем использовать неявную схему: в пределах всего шага по времени (¿0 < t < ¿0 + Д^) значения Ф и р принимаются равными новым значениям.

Тогда уравнение (12) аппроксимируется следующим дискретным аналогом:

(р °Ф° - р°РФ°) Д + ( Зп - З-) АфАг + (Л, - Зе) АгАг + (Л - Зл) Дг Аф = ФА У, (13)

где Ф — среднее по контрольному объему (новое, в момент времени ¿0 + АЪ) значение Ф, А V = —-АгАфАг. Отметим, что при гР = —- выполнено АУ = гРАгАфАг. Аналогично интегрируем уравнение неразрывности (8):

ш

V

to+At

I т(гр)с

to+At

|

¿Г Сф ¿X+ д

/ Ш {дгт (Григ) + ^иф) + д^ 0 СгСфСг

¿о V

Введем обозначения

(14)

ССЪ = 0.

Рг = триг, Рф = риф, Рг = григ.

(15)

При этом Рп и — значения Рг на гранях п и в, Рш и — значения Рф на гранях т и е , и Р^ — значения Рг на гранях и и С.

Учитывая (15), составим дискретный аналог уравнения (14):

(р р - рр) Ду + ( Рп - Р-) АфАг + (Рш - Ре) АгАх + (Ри - Ра) АгАф = 0. (16) Умножим (16) на ФР и вычтем из (13). Тогда

(Фр - Ф°р)

рр ау ~аГ

+ ( Зп - РпФ°) АфА-г - (Л - Р8ФР) ДфДг +

+ ( Зш - РшФ°) ДгДг - (Зе - РеФ°) ДгДг + (Л - РиФ°) АгАф --( Зл - Р,Ф°) АгАф = ФАV.

(17)

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

3. Точные решения уравнений сохранения

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

I. Рассмотрим установившееся течение вдоль оси аппликат:

д < ^ д / дФ\ дд~г{-Три Ф - д~АГГ дф) =°,

д

— (гри,) =°.

(18)

Решаем систему уравнений (18):

дФ д /дФ\

Гр Uz---г— г —

д д

д / ч

г Ж = 0-

или

дФ д дФ PUz— - — (г — ) =0,

д - д

д

puz = const.

В предположении, что Ф = Ф(^), Ф^. = Ф'ф = 0 и Г = const, получаем краевую задачу для обыкновенного дифференциального уравнения в области 0 < z < L :

// Puz ф,

Ф'' _

zz г

0, Ф(0) = Ф0, Ф(Ь) = Ф L.

Решая задачу (19), найдем:

Ф = Сгеpuzz/г + С2,

iCi + С2 = Фо, 1 Ci е pu*L/r + С2 = Фь.

(19)

Обозначая Р7

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

P UzL Г

(число Пекле), получим

Ci

Ф L - Фо

е - 1

С2 = Ф0

Ф L - Фо

е р* - 1

Фое- Ф l ер* - 1 ■

Отсюда

Фь - Фо ^pzz/l + Фоер - Фь

Ф

или

е - 1

Ф - Фо Ф L - Фо

+

е - 1

,PZz/L _ 1

е р* - 1 '

(20)

II. Рассмотрим установившееся течение в направлении полярного угла (вдоль дуги L0L при г = const, рис. 3):

А( фа А(ГдФ\

дф ф дф V г дф /

0,

(21)

0.

Считаем, что искомое решение ищется при незначительных изменениях координаты г. Поэтому полагаем г = const. Данное ограничение вполне допустимо, так как найденные нами решения будут использованы в пределах малых ячеек расчетной области.

0

Рис. 3. Область решения одномерной задачи в направлении ф

Решаем систему уравнений (21):

дФ д /Г дФ\ _

Р ф дф дф V г дф ) '

ри„

const.

В предположении, что Ф = Ф(ф), Ф'г = Ф^ = 0 и Г = const, получаем краевую задачу для обыкновенного дифференциального уравнения в области ф0 < ф < ф1 (в пределах изменения параметра дуги I = г • ф : L0 <1<L):

ри

Ффф -г^Фф = 0, Ф(фо) = Фо, Ф(ф1) = Ф1.

Г

(22)

Решая задачу (22), найдем:

Ф = С1 егрМфф/г + с2,

(Cie гр«ф фо/г + с<2 = Фо,

]с1е ^ф1/г + с2 = Ф1.

Если принять ф0 = 0, то система уравнений относительно C1 и C2 упрощается:

{

C1 + C2 = Фо'

C1 е ^ф1/г + C2 = Ф1.

ри

Обозначая Рф = г ——ф1 (число Пекле), получим

Г

C1

Ф1 - Фо

е — 1

C2 = Фо -

Ф1 - Фо Фо- Ф е — 1 е — 1

Отсюда

или

Ф = Ф1 - Фо ерффм + Фое^ - Ф1

Ф Pi + Р т '

е ^ - 1 е ^^ - 1

Ф - Фо е ^Фф/ф1 - 1

:>Рт

Ф1 - Фо е- 1 Отметим, что решение (23) будет верным при условии г = const.

(23)

1

III. Рассмотрим установившееся течение вдоль полярной оси

д , д

д(rpurФ) - д (гГ^ ) = 0,

дг дг у дг J (24)

^ (rpu) =0.

\ R„

\ R

Аг

г

Рис. 4. Область решения одномерной задачи в направлении г

Решаем систему уравнений (24):

дФ д ( ЯФ\

григ---— гГ — = 0,

д д

григ = к = const.

Рассматривая течение как одномерное, получаем уравнение относительно функции одной переменной Ф = Ф(г). В предположении Г = const задача сводится к решению краевой задачи для обыкновенного дифференциального уравнения в области R0 < г < R:

кФ'г = ГФ'г + гГФ%, Ф( Ro) = Фo, Ф( R) = Фд, где к = григ. (25)

Задача (25) приводится к виду

кл к — Г

Ф;; = - Ф;, Ф( Ro) = Фо, Ф(R) = Фд, где h = —-. г I

Общим решением дифференциального уравнения является функция Ф(г) = С1^-т-1+С2, где h + - = £ = ~ Ci С- к григ

Пусть С = fcin = "Г", к2 = г = ~г~,тогда

Ф(г) = Cirк2 + С2 = Cieк2 lnr + С2. Учитывая начальные условия, запишем систему для С и С2:

!

Ci ек2 lnRo + С2 = Фо,

к2 ln R

Ci ек2 lnR + С2 = ФR.

Находим решение задачи (25):

Фд - Фо

О,

g к2 ln R _ gfc2 ln Ro '

Отсюда

с = Ф _ ФR - Фо ек2 in Ro = Фрек2 lnR - ФReк2 lnRo

2 0 gk2 lnR _ gk2 lnRo gk2 lnR _ gk2 lnRo '

ФR - Фо 3k2 lnr , Фоеk2 lnR - ФRek2 lnRo

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

или

Ф(Г) = __^0_ек2 1ПГ- +

^ ' g к2 1пК _ gfc2 1п Ко gfc2 1пК _ gfc2 1п Ro

Ф„ _ Ф0 Ф0ек2 1п(К/Ко) _ Ф

Ф(г) = Фд Ф0 р к2 1п(г/Ко) + _^К

^ е к21п(К/Ко) — 1 ^ е к21п(К/ Ко) — 1 '

Последнее равенство можно записать в виде, аналогичном (20) и (23):

Ф — Ф0 е к2 1п(»"/ До) — 1 ФК — Ф0 = е к2 1п(К/ Ко) — 1 •

4. Физический смысл точных решений

(26)

Проще всего пояснить физический смысл полученных точных решений на примере решения (20). Исследуем зависимость Ф(г) при различных значениях числа Пекле Pz.

I. В пределе при Pz ^ 0 (в случае «чистой» диффузии) получаем

п- е P*z/L - 1 Pzz/L z ^ ^ _ ^

lim —p-— = lim ' = -, или Ф = Фр + (Ф: - Фо)-,

Pz^о ePz - 1 Pz^о Pz L L

то есть Ф( ) линейно зависит от .

II. При Pz ^ 1 для значений z, отличных от нуля, имеем

Pz / L - 1

е P - 1

_ е Pz(z-L)/L = е-Pz(L-z)/L

Ф

e-Pz(l-^l + (l - e-Pz(l-z)/l) Фо

то есть, чем больше Pz, тем больше влияние Фо на грани контрольного объема. III. При Pz ^ -1 получаем

еPzz/L - 1 e-^Pz^/L - 1 WU/L

_ = _ ~ 1 _ p-IPzlz/L

e Pz — 1 = e-|Pz I- 1 ~ 1 6

Ф= (1 - e-IPz1 z/L) ФL + e -IPzlz/LФC

то есть, при отрицательных значениях Рг чем больше Р|, тем больше влияние Ф^.

IV. Производная Ф^ (если рассматривать Ф как функцию одной переменной) имеет

вид:

- = - Фо + ePzZ/L - 1 (Ф: - Фо) ] = ФL - Фо ePz*/L • P

dz dz Фо + eP*~ 1 (Ф: Фо) eP*~ 1е L

и

и

На грани контрольного объема, расположенной посередине между узлами (то есть при г = Ь/2), имеем

= е*/2 ■ ^ Фь - Фо Р

з Р* - 1 L L ePz/2 - e-Pz/2'

Таким образом, при больших значениях IPZ| производная ФZ(L/2) близка к нулю (диффузия почти отсутствует).

Физический смысл точного решения (23) легко пояснить, исследуя аналогичным образом зависимость Ф(ф) при различных значениях числа Пекле Рф. При этом и результаты получаются точно такие же. Немногим отличается исследование зависимости Ф(г) в (26). В частности, при к2 ^ 0, что равносильно условию Pr ^ 0, получаем

е fc2 in(r/ Ro) - i lnr - 1пДо

lim

fc2^o efc2 ln(R/R) - 1 ln R — ln R

о

В этом случае зависимость Ф(г) нелинейная. Однако, как и для Ф(г), при к2 ^ 1 выполняется

ек2 1п(г/ До) _ 1

g к2 ln(R/Ro) _ 1

_ „ k2(ln(r/ R0)-ln(R/R0)) _ _-к2 ln(R/r)

о — о

Ф(г) = Ф0 + (Фд - Фо)е-к21п(Д/г).

Таким образом, на грани контрольного объема более существенно влияние Ф0. При к2 ^ — 1 имеем

Ф — Ф0 е к2 1п(Г/До) — 1

ФД — Ф0 е к2 1п(Д/ До) — 1

и это означает, что на грани контрольного объема значение Ф(г) близко к ФД.

и

5. Экспоненциальная схема

Используя полученные аналитические решения (20), (23) и (26), построим дискретный аналог для обобщенного дифференциального уравнения.

I. Найдем приближение для значения потока Р на верхней грани и контрольного объема.

В качестве профиля между узлами Р и и используем точное решение (20), заменив Ф0 на Фд, Фь на Ф^, а Р на расстояние (6г)и между точками Р и и. Тогда

Ф — ФР +

Фи - ФР

е Ри — 1

(e(P-Z)/(5z)u - 1) , Ри —

(puz)ui6z)и

Ги

0 < г < (6z)

Г 6z

Обозначим: DZ — —, mZ — —. На верхней грани контрольного объема имеем Fu

6z Z

( ч п Ги

— (rpUz)u, Du — , mu

6

Поток Ju можно представить в виде

Здесь zu — расстояние до грани и от узла Р.

J — F

ФР +

ФИ - ФР ( рРи/п е Ри — 1

Ри/ти - 1 ^

- гГ„

Фи - Ф Р Р

0Ри/т,и

Ри

- 1 6 u

где Гад = , Ри = ^^(Ь<

Следовательно,

"и 1 и

а

Фр +

Ф^ — Ф

еРи — 1

^ (еРи/ти _ 1

_ р фи — ФР рри/тп

и р

е— 1

1.

Фр +

и

Фр — Ф

е Ри — 1

(27)

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

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

Аналогично на нижней грани d контрольного объема выполняется равенство

= 1

Ф в +

р

Фр — Ф е р*~ 1

(28)

{Ьх)а, (Ьг— расстояние между узлами Р и В.

Г я

где Рл = (гри*)^ Рл

II. Найдем приближение для потока на «западной» грани т контрольного объема.

В качестве профиля между узлами Р и Ш используем «средневзвешенное» по этой грани значение Ф на основе решения (23), заменив Ф0 на Фр, Ф1 на Ф^, а ф1 на (Ьф)ад (см. рис. 2). Тогда

Ф = Фр +

Ф^ — Ф

е— 1

р ф)/(6Ф)» — ^ ,

Р'Ш = Г-Г--, 0 < ф < (Ьф^

Г т

Обозначим: В

Г Ьф

ф = -—, тф = —. На «западной» грани контрольного объема = гЬф ф ^ ^

п о р

= ———;—, = --. Здесь — угол до грани т от узла Р.

Поток можно представить в виде

/ = Р

1 ад

Фр +

Фш - Ф

_ 1

р ^рп, /тт _ ^

Г» Фw — Фр Р'Ш

г е р™ — 1 (Ьф)

где = гВш (Ьф)ш, Рад = т^ (Ьф)ш = И^.

Г ад Вад

Следовательно,

/ = Р

•-'ад 1 ад

Ф^ — Ф; е — 1

Фр + ^ — Ф/ — 1

_ Р ФШ — ФР рад /Шад =

1 ад п с

е— 1

Р ад

Фр +

Ф Р - Ф

Р — Ф\¥ — 1

(29)

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

На «восточной» грани контрольного объема выполняется:

Фр — Фр

Ф Е +

е 1

(30)

где Fe = (риф) . Pe = r

{риф^ Ге

(бф) , (бф) — угол между узлами P и Е.

III. Найдем приближения для потока Jr на «северной» грани п контрольного объема. В качестве профиля между узлами P и N используем точное решение (26), заменив Фо на ФР, на , До на RP, R на Rn, а Ar на (бг) = Rn - RP (см. рис. 2). При этом RP < г < Rn = RP + (бг) . Тогда

Ф = ФP +

Ф n - Ф

P

еPn — 1

^e(Pn/L) ln(г/RP) - ^

ln (1 + (бг) jRp) = • L

\ / In

Г ! Км р (ри0,

где Ь = 1п —-, Рп = гп '

Г п

Г 1п(1 + Ьг / КР) с

Обозначим: Вг = -—---————, тг = —-—-————, Ьг — расстояние меж-

г 1п(1 + Ь г/КР) г 1п( г/КР) Р

ду соседними узлами вдоль полярной оси. На «северной» грани контрольного объема

Г Ь

величины Ег, Ог и тг принимают значения Еп = гп(риг) , Оп = —, тп = -——---.

п Ь 1п(гп/КР)

Здесь гп — значение полярного радиуса, соответствующее грани п. Поток /п можно представить в виде

J = F

n n

Фр + *»Z*P - 1)

е Ри — 1 V /

_г Г ^ - ФР Pn 1 Рп/тп

' n1- n р Т ° )

е Рп — 1 L г„

где Гп = DnL, Pn

n n

(pur)

Г

n

Fn

nL = — Dn

Следовательно,

J = F

n n

ф р + (?р-/-~ -1)

р рп _ 1 V /

F

n

^ - ФР еРп/шп

еРп — 1

F

n

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

Фр +

Ф P - Ф

n

зр™ - 1

(31)

Как и в предыдущих случаях, не зависит от расположения границы раздела между соседними узлами (тп в выражение для потока /п не входит).

Аналогично на «южной» грани контрольного объема выполняется:

Jя Fs

Ф S +

Ф S - Ф

P

3 р°- 1

(32)

где Fs = Гя(риг)я, Pя = Гя Г

Г

P ив.

(риг),

ln + (б г) я/Rsj, (бг) я — расстояние между узлами

Подставим выражения (27)-(32) для суммарных потоков на гранях контрольного

объема в уравнение (17).

(Фр — ФД)

рД А V

дг

+ к

Фд +

ф р - ф

N

еРп — 1

- К

- К!

Фс + Фр +

Ф с Ф

р

е р — 1 ФЕ — Фр е — 1

К*

— К,Фр АфАг + К

— АгАг + ( Ки

Фр — Фр

Фд +

Фд +

— К„Фр^ АфАг — Фр — Ф^

С>Ри,

1

и

Фр — Ф

е Ри — 1

— Еш Фр^ АгАх —

— КиФР ) АгАф —

Фр +

е р; - 1

О

— К*ФР АгАф = ФА V.

Среднее значение источникового члена Ф представим в линеаризованном виде: Ф = = Фс + ФрФр, где коэффициенты Фс и Фр определяются видом дифференциального уравнения.

Приводя подобные члены, имеем:

/рД А V 1 К; V А* '

Л

+

е Рп — 1

АфАг +

е р* - 1

К, АфАг +

е р™ — 1

АгАг +

+

рД А V А*

е Ре

р Ке АгАг + р еРе — 1 е Ри — 1

А гАф +

Л

Р; _ 1

К* АгАф — ФР АV ФР —

Ф р —

К

;

еРп — 1

АфАг Ф^у--р--К, АфАг Фс--р ™ АгАг Ф^ —

е - 1

Ке АгАх ФЕ —

К,

е Ри — 1

е р — 1 АгАф Фи —

е р; - 1

е — 1

К* АгАф Фр = Фс АУ.

Получаем следующий вид дискретного аналога:

ар Фр = аN Ф N + ас Фс + а№ Ф^ + ар Фе + аи Фи + ар Фр + Ь,

(33)

где

ар = ар + аN + ас + а№ + аЕ + аи + ар — Фр АV,

К

;

,р8

аЕ

р

ер" — 1

АфАг, ас

ре

Ке АгАг, аи =

е р — 1

Ки

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

К, АфАг, а^

АгАф, ар

е р™ — 1

А А ,

р;

е ре — 1 ' ер1 — 1 ' ер; — 1

К* АгАф,

рр А V

А

, Ь = аир Фр + Фс А V.

Значения чисел Пекле Р;, Р,, Рт, Ре, Ри, Р*, «расходы» К,, , Ке, Ки, К;* и

«проводимости» О,, , Ое, Ои, О* вычисляются по формулам, полученным в ходе построения дискретного аналога. Приведенную схему, использующую точные (экспоненциальные) решения, называют экспоненциальной.

6. Общая форма дискретного аналога

Для дискретного аналога (33) можно получить более общую формулировку, если воспользоваться исследованием С. Патанкара [5, с. 78]. Представим, например, суммарный поток Зга как среднее взвешенное значений обобщенной переменной в соседних

узлах:

или

где

Jn

Fn +

F

-L- Г).

ePn — 1

Ф P -

F

-L- n

ePn — 1

Ф

n

Jn

j* = Dt = b фр - a фn ,

Dn

А

P

n

=р1

В = Pn +

P

• n

(34)

e Pn — 1

А + Pn.

При изменении направления координатной оси на противоположное меняется знак Рп, а коэффициенты А и В меняются местами. Тогда

- J* = В * Ф]у - А* Ф

P,

А* = А(—Рп) = В(Рп), В* = В (—Рп) = А(Рп).

Таким образом, при противоположных направлениях координатной оси, то есть при значениях Рп разных знаков, выражения для коэффициентов А и В различаются. Решения (34) получены при условии, что направление вектора скорости и совпадает с направлением полярной оси (при Рп > 0). Если и и полярная ось противоположно направлены, то Р„ < 0. Тогда

А(Рга) = А(—|РП|) = В(|Р„|) = А(|Рга|) + |Рга| = А(|Рга|) — Рп.

Соответственно, при Рп > 0 выполнено равенство А(Рп) — А(|Р^|). Введем оператор

А(Р) = А(|Р |) +шах(—Р ;0),

| Р|

где A(|P|)

При этом

е1р I — 1

В (Р) = А(|Р |) +шах(Р ;0). Коэффициент ам запишем, используя (35):

(35)

(36)

F

n

ün =

е Рп — 1

АфА z = DnA(Pn)AфAz= Dn A(|Pn|) + max(-Fn;0) АфАг

Аналогично получаем выражение для аз, используя (36):

S

е р'~ 1

Fs АфАг = DsB(Pя)AфAz= Ds A(|Ps|) + max(Fs; 0) АфАг.

Применение операторов А( Р) и В( Р) позволяет получить обобщенную форму записи и для других коэффициентов. Тогда окончательно дискретный аналог принимает вид:

ар Фр = ам Ф^ + а3 Ф^ + а№ Ф\¥ + аЕ Фе + аи Фи + а0 Фс + Ь, (37)

где

аР — а°Р + a,N + as + aw + aE + аи + aD - ФР AV,

aw — DnA(Pn)AфAz, as — Ds В(Ps) AyAz, aw — DWA(PW)ArAz,

aE — De В(Pe)ArAz, «и — Du A(Pu) ^ ArAф, aD — DdB(Pd) rd ArAф,

a°P — ^Р^^, b — a°p ФР + AV,

а операторы A(P) и В(P) задаются формулами (35, 36).

Полученный дискретный аналог (37) обобщенного дифференциального уравнения (5), записанного в цилиндрической системе координат (9), может быть использован для численного моделирования конвективного тепломассопереноса в условиях ламинарных и турбулентных течений вязкой несжимаемой жидкости в замкнутом объеме. Вопросы реализации данной модели, касающиеся расчета поля скоростей, использования граничных условий различного рода, остались за рамками данной статьи, поскольку здесь решалась задача получения экспоненциальной схемы применительно к цилиндрическим координатам. Учет особенностей реализации модели необходим для корректной аппроксимации функции Ф между узловыми точками, для определения вида операторов A( P) и В( P) и формул расчета чисел Пекле для переменных г, ф иг. Использование дискретных аналогов, построенных на основе точных решений уравнений переноса, может оказаться более предпочтительным с учетом современного уровня вычислительной техники.

СПИСОК ЛИТЕРАТУРЫ

1. Белова, О. В. Метод контрольного объема для расчета гидравлических сетей / О. В. Белова, В. Ю. Волков, А. П. Скибин // Инженерный журнал: наука и инновации. — 2013. — Вып. 5. — C. 1-14.

2. Берковский, Б. М. Вычислительный эксперимент в конвекции / Б. М. Берковский, В. К. Полевиков. — Минск : Университетское, 1988. — 167 с.

3. Будак, Б. М. Кратные интегралы и ряды / Б. М. Будак, С. В. Фомин. — М. : Наука, 1965. — 608 с.

4. Патанкар, С. В. Численное решение задач теплопроводности и конвективного теплообмена при течении в каналах / С. В. Патанкар. — М. : Изд-во МЭИ, 2003. — 312 с.

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

6. Флетчер, К. Вычислительные методы в динамике жидкостей / К. Флетчер. — М. : Мир, 1991. — Т. 1. — 502 с.

7. Heiss, A. Numeris^ und experimentelle и^егБисЬип£еп der laminaren und turbulenten Konvektion in einem ges^lossenen Behälter. Dissertation / A. Heiss. — Münzen : Lehrstuhl A fur Thermodynamik, Te^nis^e Universität Munden, 1987. — 203 p.

REFERENCES

1. Bеlova O.V., Volkov V.Yu., Skibin A.P. Mеtod kontrolnogo obyеma dlya raschеta gidravlichеskikh sеtеy [Control Volume Method for Cakulation of Hydra^te Networks]. Inzhsnsmyy zhurnal: nauka i innovatsii, 2013, iss. 5, pp. 1-14.

2. Bеrkovskiy B.M., Potevikov V.K. Vychislitßlnyy ekspßrimßnt v konvßktsii [Computing Experiment in Convedtion]. Minsk, Univеrsitеtskoе Publ., 1988. 167 p.

3. Budak B.M., Fomin S.V. Kratnyß intßgraly i ryady [Multiple Integrals and Series]. Mos^w, Nauka Publ., 1965. 608 p.

4. Patankar S.V. Chislennoe reshenie zadach teploprovodnosti i konvektivnogo teploobmena pri techenii v kanalakh [Computation of Conduction and Duct Flow Heat Transfer]. Moscow, Izd-vo MEI Publ., 2003. 312 p.

5. Patankar S. Chislennye metody resheniya zadach teploobmena i dinamiki zhidkosti [Numerical Heat Transfer and Fluid Flow]. Moscow, Energoatomizdat Publ., 1984. 152 p.

6. Fletcher K. Vychislitelnye metody v dinamike zhidkostey [Computational Methods for Fluid Dynamics]. Moscow, Mir Publ., 1991, vol. 1. 502 p.

7. Heiss A. Numerische und experimentelle Untersuchungen der laminaren und turbulenten Konvektion in einem geschlossenen Behalter. Dissertation. München, Lehrstuhl A für Thermodynamik, Technische Universitat Munchen, 1987. 203 p.

DISCRETIZATION OF DIFFERENTIAL EQUATIONS OF CONVECTION AND DIFFUSION BASED ON CONTROL VOLUME METHOD

Alexander Viktorovich Bokov

Candidate of Physical and Mathematical Sciences, Associate Professor, Department of Natural Science and Humanitarian Disciplines, Pyatigorsk Branch of Plekhanov Russian University of Economics av_bokov@mail.ru

Kuchury St., 8, 357500 Pyatigorsk, Russian Federation

Aleksey Aleksandrovich Klyachin

Doctor of Physical and Mathematical Sciences,

Head Of Department of Mathematical Analysis and Function Theory,

Volgograd State University

klyachin-aa@yandex.ru, matf@volsu.ru

Prosp. Universitetsky, 100, 400062 Volgograd, Russian Federation

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

Marina Aleksandrovna Korytova

Candidate of Physical and Mathematical Sciences, Associate Professor, Department of Mathematical Analysis and Methods of Teaching Mathematics, South Ural State University (National Research University) sma.59@mail.ru

Prosp. Lenina, 76, 454080 Chelyabinsk, Russian Federation

Abstract. The object of research is a mathematical model of convective heat and mass transfer. This model is based on the equations of hydrodynamics. Taking into account that the equations governing heat transfer, mass transfer and fluid flow are similar in many respects, and the dependent variables of our interest satisfy a generalized conservation law, we apply numerical methods to a differential equation of a generalized variable. The research objective is construction of a discrete analogue of the generalized differential equation describing the proposed mathematical model of convection in a viscous incompressible fluid.

To implement the mathematical model, it is suggested to use the control volume method. This method has a number of significant advantages in comparison to the method of finite differences. It provides strict compliance with the laws of mass, momentum, energy conservation for any group of cells, and therefore,

in the whole computational domain. The method itself and its application are sufficiently described in the book by S. Patankar [5]. The computational domain is divided into a multiplicity of control volumes. The differential equation is integrated over the control volumes. In this case, the assumption about the form of the function describing the change in the variable between two adjacent nodes is made. To represent the profiles of functions between the nodes of the grid, we use templates that approximate the exact solution obtained analytically.

As a rule, researchers use discrete analogues obtained for rectangular grids in numerical simulation. Use of such templates for curvilinear coordinate systems leads to large computational errors. Many applied problems are easier to solve in the curvilinear coordinates, for example, in cylindrical coordinates. Further, to obtain the exact templates, we need different solutions that take into account specific features of curvilinear grids. Modern computers allow performance of complex mathematical calculations. It also makes possible the use of more accurate functional dependencies for approximations. This all led to the desire of authors to get a more efficient tool — the calculated schemes and algorithms that use an exponential scheme for the cylindrical coordinate system as a basis.

The task is to obtain the best approximations of the generalized variable profiles in a variety of directions. For this purpose, we find the exact solutions of the equations of conservation separately for each coordinate. The article investigates the meaning of exact solutions made in order to ensure their correctness.

Using the obtained analytical solutions, we have built the discrete analogue of the generalized differential equation. We obtained a more general formulation for the discrete analogue based on the research of S. Patankar. The resulting discrete analogue of the generalized differential equation written in cylindrical coordinates can be used for numerical simulations of convective heat and mass transfer under conditions of axisymmetric laminar and turbulent flows of viscous incompressible fluid and fluid flows in channels and pipes. Using discrete analogues constructed on the basis of exact transfer equations may be more advantageous taking into account the present level of computer technology.

Key words: convection, heat and mass transfer, mathematical model, generalized differential equation, discrete analogue, control volume, function approximation.

i Надоели баннеры? Вы всегда можете отключить рекламу.